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ABSTRACT 


This  report  documents  a  computer  program, 
SWATHO,  which  optimizes  SWATH  ship  geometry, 
resulting  in  a  minimum  resistance  and  EHP  for 
a  specified  speed  of  the  ship.  Contained  in 
the  report  is  a  detailed  description  of  the 
procedure  for  assembling  the  input  deck  for  use 
on  the  CDC-6000  computer  system  at  DTNSRDC,  running 
the  program  on  this  system  and  understanding  the 
output  from  the  program.  Also  included  are  main 
program  and  subroutine  descriptions  as  well  as  the 
theory  behind  the  computations  being  performed. 


ADMINISTRATIVE  INFORMATION 

This  investigation  was  authorized  under  a  direct  funded  block  from  the 
Naval  Material  Command  (NAVMAT  08T2)  under  Program  Element  62543N,  Task 
Area  ZF43-421-001,  and  administered  bv  the  David  W.  Taylor  Naval  Ship 
Research  and  Development  Center  (DTNSRDC) ,  Work  Unit  1500-102-30,  and 
1-1500-104-70. 


INTRODUCTION 

The  computer  program  SWATHO,  developed  at  DTNSRDC,  optimizes  SWATH  ship 
geometry  to  arrive  at  minimum  values  of  resistance  and  effective  horsepower  for 
a  predetermined  speed.  The  minimum  values  arrived  at  by  the  optimization  are 
a  result  of  initial  ship  geometry  and  constraints  on  the  geometry  as  specified 
by  the  program  user. 

The  calculations  use  standard  EHP  prediction  methods  accounting  for 
frictional  resistance,  form  drag,  and  wave  resistance.  Wave  making  predictions 
are  made  using  the  work  of  Lin  and  Day"''  on  twin  hull  ships. 

Contained  in  this  report  is  a  user  manual  which  explains  how  to  assemble 
the  input  deck  as  well  as  how  to  run  the  program  on  the  CDC-6000  computer 
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system  at  DTNSRDC.  A  description  of  the  output  and  how  to  interpret  it  is  also 
included  along  with  art  overview  of  how  the  main  program  works.  Following  the 
user's  manual,  a  maintenance  manual  is  presented,  which  contains  a  detailed 
descrip t i on  ot  the  common  blocks,  functions  and  subroutines,  and  a  tree 
diagram.  The  appendices  contain  a  description  of  the  objective  and  penalty  func¬ 
tions  used  for  minimization  of  effective  power,  a  presentation  of  the  theory 
behind  the  computational  procedures  for  the  resistance  predictions,  and  an  explana 
tion  of  the  relationships  between  Chebychev  series  and  SWATH  hullform  coefficients 
This  program  will  assist  naval  architects  in  finding  a  minimum  resistance 
hull  design  which  satisfies  the  given  geometric  constraints.  The  User's  Manual 
portion  of  this  document  is  intended  to  provide  the  user  with  the  information 
needed  to  prepare  input  for  the  program,  and  with  a  description  of  the  resulting 
output.  The  Maintenance  Manual  provides  the  documentation  needed  to  understand 
the  functioning  of  the  program  and  its  various  subprograms,  at  a  detailed  enough 
level  that  changes  could  be  made  in  the  program  if  necessarv. 
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USER'S  MANUAL 


PROGRAM  DESCRIPTION 

The  Program  SWATHO  optimizes  a  SWATH  ship  geometry  to  obtain  a  minimum 
value  of  resistance  and  EHP  for  a  specified  speed  of  the  ship.  The  program 
can  be  better  understood  if  it  is  divided  into  four  parts,  each  part  having 
a  specific  purpose  and  result. 

The  first  part  of  the  program  reads  the  input  data  that  the  user  has 
prepared  on  cards;  thus  initializing  the  ship  description  as  well  as  setting 
up  the  constraints.  The  array  of  scaling  factors,  which  is  used  as  the  initial 
geometry,  is  optimized  and  initially  set  to  1  in  this  part.  After  each  run, 
the  values  in  the  scaling  array  which  correspond  to  the  value  of  minimum  EHP 
at  that  point  are  input  to  the  next  run  of  the  program,  keeping  the  rest  of 
the  INPUT  cards  the  same.  The  actual  description  and  format  of  the  INPUT 
cards  are  found  in  the  INPUT  DESCRIPTION  and  INPUT  EXAMPLE  of  this  USER'S  MANUAL. 

The  second  part  of  the  program  actually  performs  the  optimization  of 
the  EHP  function.  A  detailed  description  of  the  method  used  and  the  actual 
subroutines  which  perform  the  optimization  can  be  found  in  the  MAINTENANCE  MANUAL. 

The  third  part  of  the  program  evaluates  the  objective  function,  thus 
finding  intermediate  values  and  a  final  value  of  EHP.  A  penalty  function 
is  also  calculated  in  this  part,  based  on  any  violation  of  constraints 
by  the  ship  geometry  at  that  time.  Then  more  optimization  is  performed. 

This  optimization  cycle  is  repeated  until  the  computational  time  is  exceeded 
or  sufficient  convergence  is  reached. 

The  last  part  of  the  program  outputs  the  initial  ship  description, 
intermediate  curves  and  a  design  description,  as  well  as  resistance  calculations 
over  a  range  of  speeds  and  other  results  of  calculations.  A  detailed 
description  of  the  output  is  found  in  the  OUTPUT  DESCRIPTION  and  OUTPUT 
EXAMPLE  in  the  USER'S  MANUAL. 


3 


INPUT  DESCRIPTION 


Diagrams  showing  the  orientation  of  the  geometry  represented  by  the  variables 
involved  in  this  input  description  are  found  in  Figure  I. 

INPUT  CARDS  (Total  number  of  input  cards  is  8) 

1.  The  first  card  contains  the  alphanumeric  title  of  the  model. 

Variable:  TITLE  (8) 

Format:  8A10 


1 _ 11 _ 21 _  31 _ 41 _ 51 _  51 _ 71 _ 

|TITLE  (Dl  TITLE  (2)  I  TITLE  (3)1  TITLE  (4)[TITLE  (5)[TITLE  (6)1  TITLE  (7 )  j  TITLE  (87 


2.  Data  on  the  second  card  are: 


XLSI  - 

Length  of  strut  in  ft 

HSI  - 

Draft  of  strut  in  ft 

TSMAXI  - 

Strut  thickness  in  ft  at  XLSI/2 

CWPI  - 

Waterplane  area  coefficient 

CLCFI  - 

Waterplane  moment  coefficient 

KG  - 

Height  of  CG  above  baseline 

Variable:  XLSI,  HSI,  TSMAXI,  CWPI,  CLCFI,  KG 
Format:  6F10.6 


1 _ 

|  XLSI 


11 


HSI 


21 _ 31 

TSMAXI  |  CWPI 


|  CLCFI 


KG 


J 


4 


a 


•  liii 


3.  Data  on  the  third  card  are: 


XLBI 

- 

Length  of  body  in  ft 

BDIAI 

- 

Bodv  diameter  in  ft  at  XLBI/2 

BSI 

- 

Separation  of  hull  centerlines 

CSTRTI 

- 

Dist.  of  strut  CL  fwd  of  body  CL 

CPI 

- 

Body  prismatic  coefficient 

CLCBI 

- 

Body  moment  coefficient 

Variable:  XLBI,  BDIAI ,  BSI,  CSTRTI,  CPI,  CLCBI 

Format  :  6F10.6 

1 _ 11 _ 21  31 _ 41  51 

XLBI  [BDIAI  1  BSI  j  CSTRTI  1  CPI  ~|  CLCBI  | 

4-5  Data  on  the  fourth  and  the  fifth  card  are: 

S  -  Array  of  scale  factors 
Variable:  S 
Format  :  8  F10.6 


71 _ 

I  S  (8) 


1 


11 


21 


31 


6.  Data  on  the  sixth  card  is: 

OPTS  PI)  -  Optimization  speed  in  knots 
Variable:  OPTSPD 
Format:  F10.6 


OPTSPD 


7.  Data  on  the  seventh  card  are: 


DISPMN 

- 

Minimum  ship  displacement  in  long  tons 

DFAC 

- 

Minimum  draft/body  diameter  ratio 

DRFTMX 

- 

Maximum  draft  of  ship 

WMAX 

- 

Maximum  width  of  ship 

XFWD 

- 

Minimum  distance  from  body  leading  edge 
to  strut  leading  edge 

XAFT 

- 

Minimum  distance  from  body  trailing  edge 
to  strut  trailing  edge 

Variable : 

DISPMN, 

DFAC,  DRFTMAX,  WMAX,  XFWD,  XAFT 

Format : 

6F10.4 

1 

11  21 

31  41  51 

f  DISPMN 

1  DFAC  1  DFTMAX 

I  WMAX  |  XFWD  |  XAFT 

8.  Data  on  the  eighth  card  are: 


MINGMT 

- 

Minimum 

TRANSVERSE  GM 

MINGML 

- 

Min imum 

LONGITUDINAL  GM 

LODMAX 

- 

Maximum 

body  length  over  diameter 

ratio 

LODMIN 

- 

Minimum 

body  length  over  diameter 

ratio 

TSMIN 

- 

Minimum 

strut  thickness  at  50  percent  chord 

Variable: 

Format: 

MINGMT, 

5F10.1 

MINGML, 

LODMAX,  LODMIN,  TSMIN 

1 _ 11 _ 21 _  3i _ 41 _ 

I  MINGMT  \  MINGML  1  LODMAX  1  LODMIN  |  TSMIN | 
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INPUT  EXAMPLE 


Example  of  Input  Data  (Optimization): 


I 

II 

21  31  41 

51 

61 

71 

TITLE  (6) 
_ 

TITLE  (7) 

TITLE  (8) 

SWATH  OPTIMIZATION  PROGRAM 


l 

11 

21 

31  41 

51 

XLSI 

HSI 

_  _ 

|  TSMAXI 

CWPI 

CLCFI 

f 

KG 

210. 

10.725 

6.500 

.8511 

.014286 

20. 

1 

11 

21 

31 

41 

51 

XLBI 

I  BDIAI 

]  BSI 

CSTRTI 

CPI 

j  CLCBI 

260. 

14.30 

68.50 

15. 

.849 

.01538 

1 

11 

21 

31 

41 

51 

61 

71 

S(l) 

S(2) 

-  - 

S(3) 

S  (4) 

j  S<5> 

S(6) 

S(7) 

.  _ 

S(8) 

1.0154 

.9900 

0.9900 

0.9783  -1.0 

.9032 

-1.0000 

1.005 

1 

11 

21 

31 

S(9) 

S(10) 

S(ll)  j 

S(12) 

.9900 

1.2500 

1.0000 

.5000 

1 _ 

OPTSPD 

20. 


1 

11 

21 

" 

31 

41 

51 

DISPMN 

DFAC 

DRFTMX 

WMAX 

XFWD 

XAFT 

2675.0 

1.67 

30.0 

106.0 

7.50 

30.0 

I  II  21  31  41 


MINGMT 

MINGML 

LODMAX 

LODMIN 

TSMIN 

30.00 

30.00 

20.00 

16.00 

5.0 
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OUTPUT  DESCRIPTION 


The  first  page  of  output  consists  of  an  echo  print  of  the  input 
including  the  S  array  of  scaling  factors.  Following  the  array  of  scaling 
factors,  NL  and  NF  (value  of  counters  in  the  optimization  procedure) 

P  (the  penalty  value),  PEHP  (the  function  to  be  minimized),  and  a  repeat 
of  tie  scaling  factors  (S  array)  are  printed. 

The  second  page  of  output  consists  of  an  intermediate  plot  of  a  body 
sectional  area  curve  and  strut  waterplane  outline  curve  where  the  symbols 
B  and  S  stand  for  body  and  strut  respectively.  At  the  bottom  of  the  plot, 
the  values  for  the  strut  and  body  geometry  are  printed  with  the  effective 
horsepower  required  for  the  design. 

The  third  page  of  output  lists  values  for  intermediate  design 
consisting  of  the  symmetric  and  antisymmetric  Chebychev  coefficients  for 
the  body  and  strut  and  the  wetted  surface  and  ship  geometry  for  the  optimization 
speed.  Also  a  table  of  speed-length  ratios,  wave  resistance,  frictional 
and  residual  resistance  coefficients  as  well  as  the  actual  resistances 
and  the  EHP  predicted  is  printed  for  the  range  of  speeds  from  10  to  26  knots. 

The  fourth  page  of  output  consists  of  a  tabular  listing  of  all  of 
the  scale  factors  of  the  S  array  corresponding  to  the  strut  and  body 
geometric  variables,  as  well  as  the  penalty  value  and  PEHP.  The  second 
difference  array  D(*)  is  printed  out  at  the  bottom  of  the  listing.  Also 
printed  is  NL,  NF,  PEHP  and  the  values  of  the  S  array  corresponding  to 
a  minimum  PEHP  so  far.  These  S  values  are  input  by  the  user  for  the  next 
run  of  the  program.  Often,  when  the  execution  time  limit  is  exceeded,  the 
output  will  stop  in  the  middle  of  this  tabular  listing,  in  which  case 
the  minimum  PEHP  is  found  and  the  S  array  of  scale  factors  are  picked 
out  from  the  table  and  used  as  input  for  the  subsequent  run. 

The  sequence  of  body  sectional  area  curve  and  strut  waterplane  curve, 
intermediate  design  description,  followed  by  tabular  output  is  repeated 
until  the  time  limit  is  exceeded  or  until  convergence  of  the  PEHP  function 
is  reached. 


When  the  program  has  reached  a  prescribed  degree  of  convergence, 
a  final  page  is  printed  out,  entitled  "Wave  Resistance  for  Strut  and  Body 
Geometric  Characteristics."  This  page  contains  the  optimum  strut  and  body 
geometric  characteristics,  separation  in  feet,  B/LB,  strut  centerline 
distance  from  body  centerline  (STRUT  OFFSET),  strut  offset/LE,  the  wave  drag 
for  strut  and  body  and  the  frictional  drag  of  body  and  strut.  The  coefficients 
of  resistance  as  well  as  the  minimum  EHP  calculated  at  certain  speeds  are  also 
given. 
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OUTPUT  EXAMPLE 
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JOB  CONTROL  FOR  SWATHO 


The  Fortran  code  file  of  SWATHO  is  stored  on  a  DTNSRDC  Computer  Center  disk 
pack.  To  access  the  disk  pack  and  make  use  of  the  program  the  following  directions 
should  be  followed. 

One  must  first  set  up  an  account  with  the  DTNSRDC  Computer  Center  and  acquire 
a  user  ID.  A  Computer  Center  Reference  Manual  will  provide  the  user  with  a 
further  explanation  of  the  Job  Setup  Control  cards. 


1.  To  get  the  program  from  the  disk  pack  to  a  permanent  file: 

JOB  CARD 
CHARGE  CARD 

PAUSE.  JOB  REQUIRES  DV4865. 

MOUNT, VSN=DV4865 , SN=CHRELIB . 

REQUEST, LOG, *PF. 

ATTACH, OLDPL,OPTSWATH, ID=CHRE, SN-CHRELIB . 

UPDATE (P,F) 

FTN(I=COMPILE,R=3) 

CATALOG, LOG, USER'S  NAME  OF  FILE,ID=USER  ID 

7/8/9 

6/7/8/ 9 

2.  To  run  the  program  which  is  stored  on  the  user's  permanent  file: 

JOB  CARD 
CHARGE  CARD 

ATTACH, USER'S  NAME  OF  FILE.HMJSER  ID. 

FIRST  SEVEN  LETTERS  OF  NAME  OF  FILE  (PL=50000) 

7/8/9 

INPUT  CARDS 
6/7/8/9 


The  user's  permanent  file  which  is  set  up  with  the  catalog  card  in  Part  1  will 
be  purged  after  30  days  of  inactivity. 
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MAINTENANCE  MANUAL 


This  part  of  the  report  contains  a  description  of  the  main  program, 

SWATHO  and  the  subroutines  directly  called  by  it  as  well  as  descriptions 
of  how  other  important  subroutines  fit  together.  Also  included  in  this  manual 
are  descriptions  of  the  common  blocks,  functions  and  subroutines  and 
a  functional  listing  of  the  functions  and  subroutines. 

INTRODUCTION 

The  main  program,  SWATHO,  first  calls  subroutine  INPUT.  Subroutine  INPUT 
reads  the  input  cards  prepared  by  the  user.  The  strut  geometric  characteristics, 
waterplane  area  and  moment  coefficients  as  well  as  the  height  of  the  CC  above  the 
baseline  are  read.  Next,  the  body  geometric  characteristics,  separation  of  hull 
centerlines,  placement  of  the  strut  centerline  in  relation  to  the  body  centerline 
and  the  body  prismatic  and  moment  coefficients  are  read.  Following  that  the  S 
array  of  scaling  factors,  the  optimization  speed  and  limiting  values  are  read. 

After  returning  from  INPUT  the  program  calls  PRAXIS,  the  optimization 
routine,  which  is  described  in  detail  in  the  following  pages  of  this  report. 

Upon  the  completion  of  the  optimization  by  PRAXIS,  the  subroutine  SCALE  is 
called  which  sets  up  an  array  of  current  working  values  of  X,  each  member  of  the 
array  corresponding  to  one  of  the  ship  geometric  characteristics.  The  X  array 
is  calculated  by  multiplying  S  by  XI  where  the  S's  are  the  scale  factors  and 
the  XI 's  are  the  initial  values  of  the  ship  geometry.  Both  S  and  XI  are  read 
by  subroutine  INPUT.  Initially  the  members  of  the  S  array  all  have  values  of  I 
and  then,  on  subsequent  runs,  the  S's  correspond  to  the  lowest  value  of  EHP  so 
far  encountered  by  the  user  in  the  output  from  the  previous  run  of  the  program. 

All  of  the  initial  variables  are  left  unchanged. 

Subroutine  CHEB,  called  by  the  main  program  next,  calculates  the  Chebychev 
coefficients  for  the  ship  as  well  as  the  wetted  surface  of  the  body  and  strut.* 

It  also  calls  another  subroutine,  PCHEB,  which  plots  the  strut  waterplane  and 
body  sectional  area  curves. 

LISTEX  is  then  called  and  values  cr  ship  geometry,  etc.  are  listed  to  display 
the  proximity  of  the  solution  to  the  extreme  limits. 


*  See  Appendix  C  for  a  discussion  of  the  relationship  between  the  Chebychev 
coefficients  and  the  hullform  coefficients. 
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Function  EHP  is  lastly  called  at  which  time  final  design  resistance  and  power 
is  output.  Whether  the  program  gets  this  far  is  contingent  upon  how  much  time 
the  computer  is  given  for  the  computations. 

Subroutine  PRAXIS  accomplishes  the  optimization  of  the  objective  function, 

EHP.  PRAXIS  finds  the  minimum  of  a  function  using  the  principal  axis  method  as 

2 

implemented  by  Brent.  The  arguments  of  the  subroutine  are  TO,  HO,  N,  IPRIN,  X, 

F,  and  FMIN.  PRAXIS  attempts  to  find  this  minimum  such  that  if  XO  is  the  true 
local  minimum  near  X,  then  NORM  (X-XO)-TO  +  SQUARE  ROOT  (MACHEP)  *NORM  (X),  where 
MACHEP  is  the  machine  precision,  the  smallest  number  such  that  1  +  MACHEP >  1,  and 
TO  is  the  tolerance. 

The  values  of  TO,  HO,  N  and  IPRIN  are  set  in  the  data  statement  in  SWATHO. 

The  value  of  MACHEP  should  be  2**-47  (about  7.105  E-15)  on  the  CDC  6000  series  for 
single  precision  arithmetic  or  16**-13  (about  2.23  E-16)  for  double  precision 
on  the  IBM  370  system.  HO  is  the  maximum  step  size,  and  should  be  set  to  about 
the  maximum  distance  from  the  initial  guess  to  the  minimum.  This  value  of  HO 
affects  the  initial  rate  of  convergence,  in  the  current  program  HO  *  .25.  N  is 
the  number  of  variables  upon  which  the  function  depends,  in  this  case  N  ■  12.  N  must 
be  at  least  2.  IPRIN  controls  the  printing  of  intermediate  results  of  the  optimi¬ 
zation. 

X  is  an  array  which  contains  a  guess  of  a  point  of  minimum  and  an  estimated 
point  of  minimum,  on  entry  and  return,  respectively.  F(X,N)  is  the  function  to 
be  minimized,  in  our  case  it  is  PEHP.  F  needs  to  be  a  real  function  declared 
EXTERNAL  in  the  calling  program.  FMIN  is  set  to  the  minimum  value  of  F  that  is 
found. 

The  approximating  quadratic  form  is 

Q(X')  =  F (X,N)  +  i  *  (X*  -X)  TRANSPOSE  *A*  (X*  -X). 

Here  A  is  INVERSE  (V-Transpose)  *D*  INVERSE  (V)  where  V(*,*)  is  the  matrix  of 
search  directions  and  D(*)  is  the  array  of  second  differences. 

PRAXIS  in  operation  proceeds  as  follows.  PRAXIS  calls  PRMIN  which  minimizes 
the  function  along  first  direction  V(*,J).  PRMIN  then  calls  function  PRLIN  which 
is  a  function  of  one  real  variable,  which  it  minimizes.  Function  RANDUM  is  called 
and  a  random  number  is  returned  to  PRAXIS  to  avoid  resolution  valleys.  PRMIN 
continues  to  minimize  the  PRLIN  function  in  different  directions.  PRQUAD  is  then 
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called  and  this  subroutine  tries  a  quadratic  extrapolation  in  case  minimization 
is  being  done  in  a  curved  valley.  PRAXIS  then  calls  PRFIT  to  find  principal 
values  and  directions  of  the  quadratic  forms  and  PRSORT  is  called  to  sort  the 
eigenvalues  and  eigenvectors.  PRPRIN  prints  out  the  scaling  factor  on  the  first 
page  of  the  output  as  is  described  in  the  output  section  of  the  USER’s  manual  and 
VCPRINT  prints  the  second  difference  array. 

The  optimization  procedure  is  performed  on  the  function  PEHP.  PEHP  gets  its 
value  from  the  EHP  calculated  by  the  function  by  that  name  and  a  penalty  P,  so 
that  PEHP  =  (1  +  P)*  EHP.  P  is  determined  from  an  array  calculated  in  subroutine 
CKLIM.  CKLIM  uses  the  constraints  determined  by  subroutine  INPUT  to  check  the 
current  ship  values  for  violation  of  these  constraints.  The  checking  procedure 
is  quantified  with  an  array  GG  of  violation  coefficients,  in  the  following  way. 

The  current  value  of  ship  displacement  is  multiplied  by  a  member  of  a  scaling 
array,  ALPHA,  given  in  a  data  statement  in  CKLIM.  This  results  in  a  violation 
coefficient  GG(1)  for  displacement.  Nineteen  of  these  coefficients  are  calculated 
and  this  array  of  GG(I)'s  is  used  in  function  PEHP.  Since  GG(I)>0  means  that 
the  constraint  is  violated  whereas  GG(I)<0  means  that  the  constraint  is  satisfied, 
the  function  PEHP  totals  up  all  of  the  GG(I)'s  greater  than  zero  and  this  becomes 
the  value  of  P,  used  to  calculate  PEHP.  For  further  details  on  penalty  functions 
see  Appendix  A. 

Other  important  subroutines  have  not  been  mentioned  yet.  Subroutine  RWAVE, 
called  by  Function  EHP,  computes  the  auxiliary  wave  resistance  functions  T  and  W 
as  explained  in  Appendix  B,  using  the  Chebychev  expansions  calculated  in  subroutine 
CHEB,  as  discussed  in  Appendix  C.  Subroutine  RINTEG  integrates  these  auxiliary 
functions  with  subroutine  INIT  initializing  the  arrays.  The  symmetric  matrices 
of  T  and  W  are  reflected  by  subroutine  REFLKT.  Using  the  Chebychev  coefficients, 
subroutine  WSURFB  determines  the  wetted  surface  of  the  body  of  revolution  and 
WSURFS  determines  the  wetted  surface  of  the  strut. 

Subroutine  ROUT  is  called  by  function  EHP,  it  computes  values  of  Froude  number 
and  speed-length  ratio  as  well  as  wave  resistance  coefficients  and  frictional  drag 
coefficients  and  returns  the  value  of  EHP  to  function  EHP.  It  also  prints  the 
intermediate  and  final  design  values.  Which  one  of  these  is  printed  is  determined 
by  the  value  of  KEY.  If  KEY  =  0  no  output  is  printed,  instead  only  the  compu¬ 
tations  in  ROUT  are  performed.  If  KEY  =  1  the  optimization  is  complete  and  just 
the  final  design  is  printed.  If  KEY  =  2  when  ROUT  gets  to  the  output  section. 


the  intermediate  results  are  printed.  The  different  values  of  KEY  are  set  in 
the  subroutines  which  call  FUNCTION  EHP  which  in  turn  calls  ROUT;  KEY  is  set  to 
0  in  PEHP,  1  in  SWATHO  and  2  in  AUXPLOT. 


Figure  i.  -  TREE  DIAGRAM  OF  SWATH  PROGRAM 
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Subroutine  Name 
AUXPLT 

BESS  J 
CHEB 

CKLIM 

INPUT 

LISTEX 

MAPRNT 

PCHEB 

PLOT  I 

PLOT2 

PLOT3 

PL0T4 

POMIT1 

POMIT2 

PRAXIS 

PRFIT 

PRMIN 

PRPRIN 


FUNCTIONAL  LISTING  OF  FUNCTIONS  AND  SUBROUTINES 


Description 

Prints  titles  and  labels  for  body  sectional 
area  and  waterplane  curves 

Evaluates  Bessel  function 

Determines  Chebychev  coefficients  for 
body  and  strut 

Checks  for  violation  of  geometric  constraints 

Reads  in  initial  values  of  ship  geometry 
and  constraints 

Lists  various  values  of  the  design  to  see  how 
close  to  the  limits  the  solution  lies 

Prints  columns  of  matrix 

Plots  body  sectional  area  curve  and  water- 
plane  outline  curve 

Sets  up  spacing  and  determines  values  of  the 
axes  for  printer  plot 

Establishes  formula  for  computing  location 
in  the  image  region  where  the  point  will  be 
plotted,  in  the  printer  plot 

Assigns  an  alpha  character  to  each  point 
to  be  plotted  on  the  printer  plot 

Prints  image  of  completed  graph  on  the 
printer 

Causes  certain  grid  lines  on  printer  plot 
to  be  deleted 

Causes  values  at  grid  lines  to  be  deleted 

Finds  the  minimum  of  the  objective  function 
using  the  principal  axis  method 

Fits  points  to  a  quadratic  curve 

Minimizes  the  objective  function,  F,  in  one 
dimension 

Prints  intermediate  results  and  calls  the 
routine  which  does  the  printer  plotting 
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FUNCTIONAL  LISTING  OF  FUNCTIONS  AND  SUBROUTINES  (CONT) 


Subroutine  Name 
PRQUAD 

PRSORT 

QPLOTZ5 

REFLECT 

RINIT 

RINTEG 

ROUT 

RWAVE 

SCALE 

SIMPSN 

SUM 

VCPRINT 

WSURFB 

WSURFS 


Description 

Looks  for  minimum  of  F  along  a  quadratic 
curve 

Sorts  the  elements  of  a  matrix  to  find  the 
eigenvalues  and  eigenvectors 

Calculates  scaling  information  to  label  the 
plot 

Reflects  matrices  of  auxiliary  wave  resistance 
functions 

Initializes  the  auxiliary  wave  resistance 
funciton  matrices  to  zero 

Evaluates  the  integrand  of  the  auxiliary 
wave  resistance  function 

Calculates  EHP,  wave  resistance 
coefficients  and  constants,  and  also 
prints  the  intermediate  and  final  design 
results 

Computes  the  auxiliary  wave  resistance 
functions 

Calculates  working  values  of  ship  geometry 

Sets  up  Simpson's  multipliers  for  numerical 
integration  by  Simpson’s  rule 

Computes  the  components  of  wave  resistance 
as  a  matrix  multiplication  involving  the 
Chebychev  coefficients  and  the  auxiliary 
wave  resistance  function  matrices 

Depending  on  the  option;  prints  the 
second  difference  array,  scale  factors, 
value  of  quadratic,  and  value  of  the 
objective  function  as  part  of  the 
optimization 

Determines  surface  area  of  a  body  of 
revolution  given  by  Chebychev  coefficients 

Determines  surface  area  of  strut  whose 
thickness  distribution  is  given  by 
Chebychev  coefficients 


Function  Name 
CFITTC 

EHP 

EVAL 

FORMDR 

PEHP 

PRLIN 

RANDUM 

YINTP 


FUNCTIONAL  LISTING  OF  FUNCTIONS  AND  SUBROUTINES  (CONT) 


Description 

Determines  frictional  resistance 
coefficient  from  ITTC  line 

Calculates  body  constants  and  frictional 
drag  coefficients 

Evaluates  the  Chebychev  series 

Evaluates  form  drag  coefficients 

The  objective  function  to  be  minimized 
by  PRAXIS 

Function  of  one  variable  minimized  by 
PRAXIS 

Returns  a  random  number  which  creates 
random  steps  during  the  optimization 

Interpolates  through  a  set  of  discrete 
data 
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COMMON  BLOCK  DESCRIPTIONS 
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COMMON /COEFS 


PURPOSE: 

NAME 

ASM 

BSM 

ABM 

BBM 

MMAX 


[ 

| 


COEFS 

strut 

stores  the 
and  body 

values  of  the  Chebychev  coefficients  for  the 

TYPE 

LENCTH 

DEFINITION 

R 

(3) 

Coefficients  of  Chebychev  Sine  Series  for  strut 

R 

(3) 

Coefficients  for  Chebychev  Cosine  Series 
for  strut 

R 

(3) 

Coefficients  of  Chebychev  Sine  Series  for 
body 

R 

(3) 

Coefficients  of  Chebychev  Cosine  Series  for 
body 

I 

Maximum  order  of  Chebychev  Series 

COMMON/ EXTLIM 

PURPOSE:  EXTLIM  defines  and  stores  physical  constants  and  geometric 
parameters 


NAME 

TYPE 

DEFINITION 

DISPMN 

R 

Minimum 

ship  displacement  in  long  tons 

DFAC 

R 

Minimum  draft/body  diameter  ratio 

DRFTMX 

R 

Max imum 

draft  of  ship 

WMAX 

R 

Max imum 

width  of  ship 

XFWD 

R 

Minimum 

distance  from  body  L.E.  to  strut  L. 

XAFT 

R 

Minimum 

distance  from  body  T.E.  to  strut  T. 

MINGMT 

R 

Minimum 

GMT  (transverse  GM) 

MINGML 

R 

Minimum 

GML  (longitudinal  GM) 

LODMAX 

R 

Maximum  length  over  diameter  ratio 

LODMIN 

R 

Minimum  length  over  diameter  ratio 

TSMIN 

R 

Minimum 

strut  thickness  at  50  percent  chord 

KG 

R 

Height  of  center  of  gravity  of  ship  above 
baseline 
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COMMON/ IN ITL 


PURPOSE: 

Store  initial  values  of 

ship  geometry 

NAME 

TYPE 

DEFINITION 

XLSI 

R 

Length  of  strut  in  ft. 

HSI 

R 

Draft  of  strut  in  ft. 

TSMAXI 

R 

Max.  strut  thickness  in  ft. 

CWPI 

R 

Waterplane  area  coefficient 

CLCFI 

R 

Waterplane  moment  coefficient 

CIYYI 

R 

Waterplane  inertia  coefficient 

CLCBI 

R 

Body  moment  coefficient 

XLBI 

R 

Length  of  body  in  ft. 

BDIAI 

R 

Body  diameter  in  ft.  at  XLB/2 

BSI 

R 

Separation  of  hull  centerlines 

CPI 

R 

Body  prismatic  coefficient 

CSTRTI 

R 

Dist.  of  strut  CL  fwd  of  body  CL 

VBI 

R 

Initial  displaced  volume  of  the  body  in 
cubic  ft. 

VS  I 

R 

Initial  displaced  volume  of  the 
strut  in  cubic  ft 

KB  I 

R 

Initial  height  of  center  of  buoyancy  above 
the  baseline 

BMLI 

R 

Initial  longitudinal  metacentric  height 

IYYI 

R 

Waterplane  moment  of  inertia 
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COMMON/NOCHG 

PURPOSE:  Checks  whether  certain  ship  geometry  variables  change  or  not. 

NAME  TYPE  LENGTH  DEFINITION 

LL  !  6  LL(I)  thru  LL(6)  represents  positions 

in  WORKVL  of  XLS,  BS,  HS,  XLB,  BDI A 
CSTRT.  In  the  data  statement  in 
SWATHO  they  are  assigned  values; 

LL(1)  =  1,  LL(2)  =  10,  LL(3)  =  2, 
LL(4)  =  8,  LL(5)  =  9,  LL(6)  =  1 


SS  R  6 


SAME  L 


Holds  array  of  scale  factors  correspon¬ 
ding  to  the  six  ship  geometry  variables 
in  the  LL  array 

Variable  which  indicates  whether  scaling 
variables  have  changed  or  not. 
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COMMON/OMEGA 


PURPOSE:  OMEGA  stores  the  values  of  variables  and  constants  for 

evaluation  of  auxiliary  wave  resistance  functions  and  ship 
resistance  coefficients 


NAME 

TYPE 

DEFINITION 

VMFPS 

R 

Speed  of  ship  (feet  per  second) 

GAMAOS 

R 

Yos 

GAMAOB 

R 

Yob 

GOSQ 

R 

(Y  )2 
os 

HSOLS 

R 

Ratio  of  draft  to  length  of  strut 

HBOLB 

R 

Ratio  of  draft  to  length  of  body 

WETS 

R 

2 

Wetted  surface  area  of  strut  (ft  ) 

WETB 

R 

Wetted  surface  area  of  body  (ft2) 

WTSURF 

R 

2 

Total  wetted  surface  area  (ft  ) 

SEP 

R 

2b  (V  1J 
os  s 

PHIS 

R 

2(h  ) 
s 

PH  IB 

R 

2<hb) 

HYob 

RATIOL 

R 

Yob/Yos 

CFS 

R 

Frictional  drag  coefficient  of  strut 

CFB 


R 


Frictional  drag  coefficient  of  body 


COMMON/PHYSCO 


PURPOSE:  PHYSCO  stores  physical  constants 


NAME 

RHO 

GNU 

G 

PI 

DELCF 


TYPE 

R 

R 

R 

R 

R 


DEFINITION 

Density  of  water 

Kinematic  viscosity  of  water 

Acceleration  due  to  gravity 

Ratio  of  circumference  to  diameter  of 
circle 

Correlation  allowance 


KTFPS  R 


Conversion  factor  for  changing  ft/sec 
into  knots 


r 

; 

i 

! 

i 


: 


COMMON /PRCOMM 


PURPOSE:  Used  by  optimization  subroutines  to  transmit  values 
amongst  themselves 


NAME 

TYPE 

DEFINITION 

FX 

R 

Function  to  be  minimized 

LDT 

R 

Used  in  step  size  computations 

DMIN 

R 

The  square  of  the  machine  precision 

NF 

1 

Function  evaluation  counter 

NT 

T 

One  dimensional  search  counter 

COMMON/PSI 

PURPOSE:  PSI  stores  the  values  of  the  variables  and  constants  used  in 
integration  routines 


NAME 

NPTSZ 

PTSAF 

EXPN 


TYPE 


DEFINITION 


I 

R 

R 


Number  of  integration  steps  from 

Y  to  y  +1 

os  os 

Scaling  factor  of  step  size  in 
integrating  from  to 

Empirical  constant  for  integration 
to  stop 


NALMAX  I 

NAL  I 

ALFA  R 


Maximum  number  of  integration  steps 

from  Y  +1  to  a 

os  max 

Counter  of  integration  steps 
Integration  variable  (a) 


ALSMAX  R 


Maximum  of  a  for  integration 


Values  of  NPTSZ,  PTSAF,  EXPN,  NALMAX  and  ALSMAX  are 
assigned  in  data  statement  in  SWATHO 


Comments: 


r 


COMMON /Q 

PURPOSE:  Used  by  optimization  subroutines  to  transmit  function 
representations  amongst  themselves 


NAME 


TYPE  DEFINITION 


V 

Q0 

Q1 


Matrix  to  define  directions  for  linear 
search 

Variables  which  define  the  plane  for  the 
quadratic  search 


QA 

QB 

QC 


R  1 

(  Variables  which  define  a  parabolic 

|  space  curve 

R  ) 


QD0 

QDl 

QF1 


l  Variables  used  in  QA,  QB,  QC  to  define 

R  j  quadratic  curve 

R  Value  of  function  defined  by  QA,  QB,  QC 
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COMMON /SPEEDS 


PURPOSE: 

NAME 

NSPEED 

SPEED 


Stores  speeds  at  which  resistance  calculations  are  made 


TYPE 

LENGTH 

DEFINITION 

I 

Number  of  speeds 

R 

20 

Speeds  in  feet  per  second 

I 


COMMON / XRPLOT  F 


PURPOSE:  XRPLOTF  stores  the  values  of  variables  for  the  plotting 
routines 


NAME 

TYPE 

DEFINITION 

XL 

R 

Value  of  abscissa  at  left-most  grid 
line 

XH 

R 

Value  of  abscissa  at  right-most  grid 
line 

YL 

R 

Value  of  ordinate  at  bottom  grid  line 

YH 

R 

Value  of  ordinate  at  top  grid  line 

XI 

R 

Increment  of  divisions  along  the 
abscissa 

YI 

R 

Increment  of  divisions  along  the 
ordinate 

YMOV 

R 

Ordinate  index  increment  number  of 
array  GRAF 

XMOV 

R 

Abscissa  index  increment  number  of 
array  GRAF 
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COMMON / XRPLOTG 


PURPOSE : 

NAME 

GRAF 


XRPLOTG  stores  the  values  and  characters  of  variables  for 
the  plotting  routines 


TYPE  LENGTH  DEFINITION 


I 


(11,204) 


Array  containing  the  image  to  be  plotted 


COMMON/  XRPLOTQ 


PURPOSE:  XRPLOTQ  stores  the  constants  and  characters  for  the  plotting 
routines 


NAME 

TYPE 

DEFINITION 

I,  II 

I 

Ordinate  scale  factor  is  10^ 

J,  JJ 

I 

Number  of  digits  following 
ordinate  decimal  point 

K,  KK 

I 

Abscissa  scale  factor  is  10 

L,  LL 

I 

Number  of  digits  following 
abscissa  decimal  point 

A,  NHL 

I 

Integer  number  of  horizontal  grid  lines 

G,  NSBH 

I 

Integer  number  of  spaces  beyond  each 
horizontal  grid  line  to  next  grid  line 

C,  NVL 

I 

Integer  number  of  vertical  grid  lines 

D,  NSBV 

I 

Integer  number  of  spaces  beyond  each 
vertical  grid  line  to  next  grid  line 

E,  HCHAR 

Horizontal  grid  character 

F,  VCHAR 

Vertical  grid  character 

M,  IX,  ISX 

I 

Number  of  horizontal  spaces 

N,  IY,  ISY 

I 

Number  of  vertical  spaces 

V 

L 

Logical  variable  =  TRUE  when  maximum 
and  minimum  values  of  the  ordinate  are 
determined 

H 

L 

Logical  variable  =  TRUE  when  maximum  and 
minimum  values  of  the  abscissa  are 
determined 
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COMMON /WORKVL 


PURPOSE: 

WORKVL  stores 

the  working  values  of  the  ship  geometry 

NAME 

TYPE 

DEFINITION 

XLS 

R 

Length  of  strut  in  ft. 

HS 

R 

Draft  of  strut  in  ft. 

TSMAX 

R 

Max  strut  thickness  in  ft. 

CWP 

R 

Waterplane  area  coefficient 

CLCF 

R 

Waterplane  moment  coefficient 

CIYY 

R 

Waterplane  inertia  coefficient 

CLCB 

R 

Body  moment  coefficient 

XLB 

R 

Length  of  body  in  ft. 

BDIA 

R 

Body  diameter  in  ft.  at  XLB/2 

BS 

R 

Separation  of  hull  centerlines 

CP 

R 

Body  prismatic  coefficient 

CSTRT 

R 

Dist.  of  strut  CL  fwd  of  body  CL 

AX 

R 

Maximum  body  cross-sectional  area  in 
sq.  ft. 

HB 

R 

Body  centerline  submergence 

DISP 

R 

Displacement  of  total  ship  in  long  tons 

AWP 

R 

Waterplane  area  of  one  strut 

NLOC 

I 

Switch  indicating  presence  of  second 
strut 

CSTRT2 

R 

Dist.  of  strut  CL  fwd  of  body  CL 

OVMFPS 

R 

Optimization  speed  in  ft/sec 
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COMMON / WORK VL  (CONT) 


NAME 

TYPE 

DEFINITION 

KB 

R 

Height  of  center  of  buoyancy  above 

baseline 

BML 

R 

Longitudinal  metacentric  height 

BMT 

R 

Transverse  metacentric  height 
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SUBROUTINE  AND  FUNCTION  DESCRIPTIONS 
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NAME: 

PURPOSE: 

CALLING  SEQUENCE: 
ARGUMENTS : 


COMMON  BLOCKS: 
SUBROUTINES  CALLED: 
CALLED  BY: 


SUBROUTINE  AUXPLT 

Prints  title  and  values  for  ship  geometry  on 
the  body  sectional  area  and  waterplane  outline 
curve.  Subroutine  PCHEB  prints  the  actual  curve 

Call  AUXPLT  (S,  N,  HP) 

S  -  Scaling  values  for  each  of  the  variables 

N  -  The  number  of  variables  upon  which  the 
function  depends 

HP  -  Horsepower 

INITL ,  NOCHG,  PHYSCO,  SPEEDS 

SCALE,  CHEB,  LISTEX,  PCHEB,  EHP 

PRPRIN 
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NAME: 


SUBROUTINE  BESSJ 


PURPOSE: 

CALLING  SEQUENCE: 
ARGUMENTS : 


COMMON  BLOCKS: 
SUBROUTINE  CALLED: 
CALLED  BY: 


SUBROUTINE  BESSJ  evaluates  the  Bessel  function 
from  order  0  to  order  N 

CALL  BESSJ  (X,  N,  VJ) 

X  -  Argument  of  the  Bessel  function 
N  -  Maximum  order  of  the  Bessel  function 
VJ  -  Array  holding  (N+l)  values  of  the  Bessel 
function  of  order  zero  up  to  N,  where 

VJ(1)  =  J  (X) 
o 


VJ (N+l)=  J  (X) 
N 

NONE 

NONE 

RINTEG,  RWAVE 
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NAME: 


FUNCTION  CFITTC 


PURPOSE : 

CALLING  SEQUENCE: 
ARGUMENT: 

COMMON  BLOCKS: 
SUBROUTINE  CALLED: 
CALLED  BY: 

COMMENTS : 


Function  CFITTC  determines  the  frictional 
resistance  coefficient  based  on  the  ITTC 
correlation  line 

C  =  CFITTC  (RN) 

RN  -  Reynolds  number  at  test  condition 

NONE 

NONE 

Function  EHP 

q  =  0  •  07  5 

(log10  (RN)  -  2J2- 


\ 

\ 


NAME: 


SUBROUTINE  CHEB 


PURPOSE: 


CALLING  SEQUENCE: 
ARGUMENTS : 

COMMON  BLOCKS: 
SUBROUTINES  CALLED: 
CALLED  BY: 


Determines  Chebychev  coefficients  for  the 
ship,  and  calculates  wetted  surface  of  body 
and  strut 

CALL  CHEB  (TITLE,  KEY) 

TITLE  -  Label  printed  on  output 

KEY  -  Variable  used  to  skip  or  include  the  printed 
output . 

COEFS,  OMEGA,  PHYSCO,  WORKVL 
WSURFB ,  WSURFS 
SWATHO,  PEHP ,  AUXPLT 
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NAME: 


SUBROUTINE  CKLIM 


PURPOSE: 


CALLING  SEQUENCE: 
ARGUMENTS : 


COMMON  BLOCKS: 


Checks  for  violation  of  constraints  and 
scales  penalties  by  the  values  in  the 
array  ALPHA 

CALL  CKLIM  (GG,  NGG) 

GG  -  actual  amount  of  violation 
NGG  -  number  of  values  to  be  checked  for 
violation 

COEFS,  EXTLIM,  WORKVL 


SUBROUTINES  CALLED: 


EVAL 


CALLED  BY: 


Function  PEHP 
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f  1 


NAME: 

PURPOSE: 


CALLING  SEQUENCE: 
ARGUMENTS : 


COMMON  BLOCKS: 


FUNCTION  EHP 

Calculates  body  constants  and  frictional 
drag  coefficients  and  calls  ROUT  to 
calculate  EHP 

EH  =  EHP  (SPEED,  TITLE,  KEY) 

SPEED  -  Speed  of  ship  in  feet  per  second 
TITLE  -  Label  printed  on  output 
KEY  -  Variable  to  skip  or  include 
the  printed  output 

NOCHG,  OMEGA,  PHYSCO,  WORKVL 


SUBROUTINES  CALLED: 


ROUT,  KWAVE,  CFITTC 


CALLED  BY: 


SWATHO,  AUXPLT ,  PEHP 


NAME: 


FUNCTION  EVAL 


PURPOSE:  Evaluates  the  Chevychev  ^ries 

MAX 

F(x)  -  E  A(M)*U(M,X)  .  B(M)*V(M,X) 

M-l 

CALLING  SEQUENCE:  EV  =  EVAL  (X,A,B,MAX) 

ARGUMENTS:  X  -  Arguments  of  the  Chebychev  series 

A  -  Coefficients  of  the  Chebychev 
Cosine  series 

B  -  Coefficients  of  the  Chevychev  Sine 
series 

MAX  -  Maximum  order  of  the  Chebychev  series 


COMMON  BLOCKS: 

NONE 

SUBROUTINES  CALLED: 

NONE 

CALLED  BY: 

CKLIM,  WSURFB,  WSURFS 

COMMENTS : 

UM(X)  -  COS  { (2M-1) (0) ’ 

VM(X)  -  SIN  2M0 

e  =  SIN-1 (X) 
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NAME: 


FUNCTION  FORMDR 


PURPOSE: 

CALLING  SEQUENCE: 
.ARGUMENT : 

COMMON  BLOCKS: 
SUBROUTINE  CALLED: 
CALLED  BY: 

COMMENTS : 


Function  FORMDR  evaluates  the  form 
drag  coefficient 

FOR  =  FORMDR  (VL) 

VL  =  Speed-length  ratio  of  ship  based 
on  strut  length 

NONE 

FUNCTION  YINTP 
ROUT 
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Figure  3.  FORM  DRAG  COEFFICIENT 


NAME: 


SUBROUTINE  INPUT 


PURPOSE: 


CALLING  SEQUENCE: 
ARGUMENTS : 


COMMON  BLOCKS: 
SUBROUTINES  CALLED: 


Reads  into  common  block/INITL/the  initial 
values  of  ship  geometry  and  constraints , 
concurrently  defining  array  XI.  It  then 
prints  the  starting  point,  along  with  an 
echo  of  the  input  values 

CALL  INPUT  (TITLE,  S,  OPTSPD) 

TITLE  -  Label  printed  on  output 

S  -  Scale  factors 

OPTSPD  -  Optimization  speed  in  knots 

EXTLIM ,  INITL 

NONE 


CALLED  BY: 


SWATHO 


NAME: 


SUBROUTINE  LISTEX 


PURPOSE: 


CALLING  SEQUENCE: 
ARGUMENTS : 

COMMON  BLOCKS: 
SUBROUTINES  CALLED: 
CALLED  BY: 


Lists  various  constraint  values  so 
that  it  can  be  determined  how  close 
to  the  extreme  limits  this  solution 
lies 

CALL  LISTEX 
NONE 

EXTLIM,  PHYSCO,  WORKVL 
NONE 

SWATHO,  AUXPLT 
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NAME: 


SUBROUTINE  MAP ENT 


PURPOSE: 

CALLING  SEQUENCE: 
ARGUMENTS : 


COMMON  BLOCKS: 


Prints  the  columns  of  the  NxN  Optimization 
matrix  V  with  a  heading  as  specified  by  OPTION 

CALL  MAPRNT  (OPTION,  V,  M,  N) 

OPTION  =  1  LABELS  NEW  DIRECTION 
2  LABELS  PRINCIPAL  AXES 

V  -  Logical  variable  =  TRUE 
M  -  The  maximum  step  size 
N  -  Number  of  variables 

NONE 


SUBROUTINES  CALLED: 


NONE 


CALLED  BY: 


PRAXIS 
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NAME: 

PURPOSE : 

CALLING  SEQUENCES: 
ARGUMENTS : 


COMMON  BLOCKS: 


SUBROUTINE  PCHEB 

Subroutine  PCHEB  plots  by  line  printer  the 
body  sectional  area  curve  and  waterplane 
outline  curve  from  the  given  Chebychev 
coefficients 

CALL  PCHEB  (AS,  BS,  AB,  BB,  NN,  TITLE) 

AS  -  Coefficients  of  Chebychev  Sine 
Series  for  strut  (Symmetric) 

BS  -  Coefficients  of  Chebychev  Cosine 
Series  for  strut  (Antisymmetric) 
AB  -  Coefficients  of  Chebychev  Sine 
Series  for  body  (Synmetric) 

BB  -  Coefficients  of  Chebychev  Cosine 

Series  for  body  (Antisymmetric) 

NN  -  Dimension  of  AS,  BS,  AB,  BB 
TITLE-  Array  containing  the  alphanumeric 

characters  of  the  title  of  the 
project. 

NONE 


SUBROUTINES  CALLED: 


PL0T1 ,  PLOT  2 ,  PL0T3 ,  PL0T4 


CALLED  BY: 


AUXPLT 
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NAME: 


FUNCTION  PEHP 


PURPOSE: 

CALLING  SEQUENCE: 
ARGUMENTS : 

COMMON  BLOCKS: 
SUBROUTINES  CALLED: 
CALLED  BY: 


PEHP  is  the  objective  function  to  be 
minimized  by  PRAXIS 

PE  =  PEHP  (S ,N) 

S  -  Array  of  scale  factors  being  optimized 
N  -  Maximum  order  of  the  Bessel  function 

OMEGA,  WORKVL 

SCALE,  CHEB ,  CKLIM ,  EHP 

PRAXIS 
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NAME : 


SUBROUTINE  PL0T1 


PURPOSE: 

CALLING  SEQUENCE: 
ARGUMENTS : 


COMMON  BLOCK: 
SUBROUTINE  CALLED: 


Subroutine  PLOT1  sets  up  spacing  and 
determines  the  values  of  the  axes 


CALL  PLOT1  (NSCALE ,  A,  B,  C,  D,  E,  F) 


NSCALE  -  Integer  array  defined  as 
follows : 


NSCALE  (1)  -  I ,  if  printed,  values  of 
the  ordinate  are  10  **  I 
times  the  actual  values 

NSCALE  (2)  -  J,  if  printed,  values  of 
the  ordinate  are  10  **  J 
times  the  actual  values 


NSCALE  (3)  -  K,  if  printed,  values  of 
the  abscissa  are  10  **  K 
times  the  actual  values 

NSCALE  (4)  -  L,  if  printed,  values  of 
the  abscissa  are  10  **  L 
times  the  actual  values 


A 

B 

C 

D 


E 

F 

XRPLOTQ 


Integer  number  of  horizontal 
grid  lines 

Integer  number  of  spaces 
between  each  horizontal 
grid  line 

Integer  number  of  vertical 
grid  lines 

Integer  number  of  spaces 
between  each  vertical  grid 
line 

Horizontal  grid  character 
Vertical  grid  character 


NONE 


CALLED  BY: 


PCHEB 
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NAME: 


SUBROUTINE  PLOT 2 


PURPOSE : 


CALLING  SEQUENCE: 
ARGUMENTS : 


COMMON  BLOCKS: 
SUBROUTINE  CALLED: 
CALLED  BY: 


Subroutine  PLOT2  examines  the  minimum  and 
maximum  values  of  the  abscissa  and  the 
ordinate  and  establishes  an  internal 
formula  for  computing  location  in  the 
image  region  corresponding  to  the  point 
to  be  plotted 

CALL  PLOT 2  (XMAX,  XMIN,  YMAX,  YMIN ,  NSCLI) 

XMAX  -  Value  of  abscissa  at  rightmost 
grid  line 

XMIN  -  Value  of  abscissa  at  leftmost 
grid  line 

YMAX  -  Value  of  ordinate  at  top  grid 
line 

YMIN  -  Value  of  ordinate  at  bottom 
grid  line 

NSCLI  -  Logical  flag  (should  be  FALSE, 
if  PL0T1  has  not  been  called 
and  standard  grid  is  desired) 

XRPLOTF,  XRPLOTQ ,  XRPLOTG 

NONE 

PCHEB 
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NAME: 


SUBROUTINE  PLOT3 


PURPOSE : 

CALLING  SEQUENCE: 

ARGUMENTS: 


COMMON  BLOCKS: 
SUBROUTINE  CALLED: 
CALLED  BY: 


Subroutine  PLOT 3  assigns  an  alpha- 
character  to  each  point  to  be  plotted 

CALL  PLOT 3  (PCHAR,  X,  Y,  SDATA,  FDATA , 
DDATA) 

PCHAR  -  Plotting  character 

X  -  Array  containing  the  X 

coordinates  to  be  plotted 
Y  -  Array  containing  the  Y 

coordinates  to  be  plotted 
SDATA  -  Integer  position  in  the  arrays 
of  the  first  ordered  pair  to 
be  plotted 

FDATA  -  =1  if  each  point  from  SDATA 
to  DDATA  is  to  be  plotted 
=2  if  every  other  point  is  to 
be  plotted 

=3  if  every  third  point  is  to 
be  plotted 

DDATA  -  Integer  position  in  the  array 
of  the  last  ordered  pair  to  be 
plotted 

XRPLOTF,  XRPLOTG 

NONE 

PCHEB 
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NAME: 


SUBROUTINE  PLOT4 


PURPOSE:  Subroutine  PLOT 4  prints  the  image  of  the 

completed  graph  on  the  printer,  including 
the  values  of  the  abscissa  and  the 
ordinate  at  the  grid  lines  outside  the 
bottom  and  left  edge  of  the  graph 


CALLING  SEQUENCE:  CALL  PLOT 4  (MCHAR,  NCHAR) 


ARGUMENTS : 

MCHAR 

Single  dimension  array 
containing  alpha  characters 
to  be  plotted  at  the  left 
of  the  graph 

NCHAR 

-  Number  of  valid  characters  in 
MCHAR 

COMMON  BLOCKS: 

XRPLOTF, 

XRPLOTG ,  XRPLOTQ 

SUBROUTINE  CALLED: 

QPL0TZ5 

CALLED  BY: 

PCHEB 
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NAME : 


PURPOSE: 


CALLING  SEQUENCE: 
ARGUMENTS : 


COMMON  BLOCKS: 
SUBROUTINES  CALLED: 


SUBROUTINE  POMIT1 

To  cause  certain  grid  .lines  on  graph  to  be 
deleted,  depending  upon  which  arguments 
are  labeled  TRUE 

CALL  P0MIT1  (T,  B,  L,  R) 

T  =  top  line  on  graph 
B  =  bottom  line  on  graph 
L  =  left  line  on  graph 
R  =  right  line  on  graph 

XRPLOTG ,  XRPLOTQ 

NONE 


CALLED  BY :  Main  program 

COMMENTS:  This  subroutine  is  part  of  a  plotting 

package  and  is  included  in  this 
program  only  to  keep  this  package 
intact 
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NAME: 


SUBROUTINE  POMIT2 


PURPOSE 

CALLING  SEQUENCE: 
ARGUMENTS : 

COMMON  BLOCKS: 
SUBROUTINES  CALLED: 
CALLED  BY: 

COMMENTS : 


Causes  horizontal  and/or  vertical  values 
at  grid  lines  to  be  deleted 

CALL  P0MIT2  (ARG) 

ARG  ARG  =  1  values  of  abscissa  at 

grid  lines  deleted 
ARG  =  2  values  of  ordinate  at 
lines  deleted 

ARG  =  3  both  sets  of  values 
deleted . 

XRPLOTQ 

NONE 

Main  program 

This  subroutine  is  part  of  a  plotting 
package  and  is  included  in  this  program 
only  to  keep  this  package  intact 
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NAME: 


SUBROUTINE  PRAXIS 


PURPOSE: 


CALLING  SEQUENCE: 
ARGUMENTS : 


COMMON  BLOCKS: 
SUBROUTINES  CALLED: 

CALLED  BY: 

COMMENTS : 


Finds  the  minimum  of  the  function 
F(X,N)  of  n  variables  using  the  principal 
axis  method.  The  gradient  of  the  function 
is  not  required. 

CALL  PRAXIS  (TO,  HO,  N,  IPRIN,  X,  F,  FMIN) 


TO 

HO 

N 


IPRIN 

X 


F 

FMIN 


is  a  tolerance 

is  the  maximum  step 

size 

the  number  of  variables 
upon  which  the  function 
depends  (must  be  at 
least  two) 

controls  the  printing 
of  intermediate 
results 

contains  oh  entry  a  guess 
of  the  point  of  minimum, 
on  return  the  estimated 
point  of  minimum 
the  function  to  be 
minimized 

is  set  to  the  minimum 
found 


PRCOMM ,  Q 


PRPRIN ,  PRMIN ,  PRQUAD ,  MAPRNT,  PRFIT , 
PRSORT ,  VCPRNT,  RANDUM,  PEHP 

SWATHO 


The  approximating  quadratic  form  is: 


s(X')  «=  F(X,N)  +  *5  *((X’-X)  -  Transpose  *A*  (X-X')) 


X  is  the  best  estimate  of  the  minimum 
A  is  Inverse  (V-Transpose)  *D*  Inverse  (V) 
V(**)  is  matrix  of  search  directions 
D(*)  is  array  of  second  differences 
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NAME: 


SUBROUTINE  PRFIT 


PURPOSE:  An  improved  version  of  MINFIT  restricted 

to  M  *  N,P  ■  0.  The  singular  values 
of  the  array  AB  are  returned  to  Q,  and 
AB  is  overwritten  with  the  orthogonal 
matrix  V  such  that  U  DIAG  (Q)  ■  AB.V, 
where  U  is  another  orthogonal  matrix. 

CALLING  SEQUENCE:  CALL  PRFIT  (M,  N,  MACHEP,  TOL,  AB,  Q) 

ARGUMENTS:  M  -  is  the  maximum  step  size 

N  -  the  number  of  variables  upon  which 
the  function  depends 

MACHEP  -  is  the  machine  precision,  the  smallest 
number  such  that  1+MACHEP>1.  MACHEP 
should  be  2**-47  (about  7.105  E-15) 
for  single  precision  arithmetic  on 
the  CDC  6000  system  or  16**-13  (about 
2.23  D-16)  for  double  precision  on 
the  IBM  370  system. 

TOL  -  Tolerance 

AB  -  Array  for  which  singular  values  are 
to  be  determined 

Q  -  Where  the  original  values  of  the  array 
are  to  be  stored. 


COMMON  BLOCKS: 

NONE 

SUBROUTINES  CALLED: 

NONE 

CALLED  BY: 

PRAXIS 
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NAME: 
PURPOSE : 


CALLING  SEQUENCE: 
ARGUMENTS : 


COMMON  BLOCKS: 
SUBROUTINES  CALLED: 
CALLED  BY: 


FUNCTION  PRLIN 

PRLIN  is  the  function  of  one  real 
variable  L  that  is  minimized  by  the 
subroutine  PRMIN 

PR  =  PRLIN  (N,  J,  L,  F,  X,  NF) 

N  -  Number  of  variable 
J  -  Defines  direction  for  search/ 
minimization  in  the  V  matrix 
L  -  The  single  variable  in  the  rotated 
coordinate  system  upon  which  the 
function  is  assumed  to  vary 
quadrat ically 

F  -  The  function  to  be  minimized 

X  -  Coordinate  of  the  initial  point 
in  the  quadratic  space 


NF  -  The  function  evaluation  counter 


Q 

NONE 

PRMIN 
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NAME: 


SUBROUTINE  PRMIN 


PURPOSE:  Minimizes  F  from  X  in  the  direction  V(*,J) 

unless  J  is  less  than  2,  when  a  quadratic 
search  is  made  in  the  plane,  defined  by 
Q0»  Ql,  X 

CALLING  SEQUENCE:  CALL  PRMIN  (N,  J,  NITS,  D2,  XI,  FI,  FK, 

F,  X,  T,  MACHEP,  H) 

ARGUMENTS:  N  -  is  the  number  of  variables  upon  which 

the  function  depends. 

J  -  Defines  direction  for  search/minimiza¬ 
tion  in  the  V  matrix. 

NITS  -  Controls  the  number  of  times  an 

attempt  will  be  made  to  halve  the 
interval 

D2  -  is  either  zero  or  an  approximation 
to  half  F 

XI  -  is  estimate  of  the  distance  from 
X  to  the  minimum  along  V  (*,  J) 

FI  -  PRLIN  (N,  J,  XI,  F,  X,  NF) 

FK  -  Logical  operator 
F  -  The  function  to  be  minimized 
X  -  Coordinate  of  the  initial  point 
in  the  quadratic  space 
T  -  Tolerance 

MACHEP  -  Machine  precision 
H  -  Step  size 

COMMON  BLOCKS:  PRCOMM,  Q 

SUBROUTINES  CALLED:  PRLIN 

CALLED  BY:  PRAXIS,  PRQUAD 
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NAME: 

SUBROUTINE  PRPRIN 

PURPOSE: 

PRPRIN  prints  intermediate  optimization 
results  including  the  number  of  function 
evaluations,  function  value,  and  current 
coordinates 

CALLING  SEQUENCE: 

CALL  PRPRIN  (N,X,  IMPRIN) 

ARGUMENTS : 

N  -  the  number  of  variables  upon  which 
the  function  depends 

X  -  contains  on  entry  a  guess  of  the 
point  of  minimum,  on  return  the 
estimated  point  of  minimum. 

IPRIN  -  controls  the  printing  of  intermediate 
results 

COMMON  BLOCKS: 

PRCOMM 

SUBROUTINES  CALLED: 

AUXPLT 

CALLED  BY: 

PRAXIS 
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NAME: 


SUBROUTINE  PRQUAD 


PURPOSE: 


CALLING  SEQUENCE: 


Looks  for  the  minimum  of  F  along  a 
quadratic  curve  defined  by  Q0,  Ql,  X, 
in  case  minimization  is  being  done  in 
a  curved  valley 

CALL  PRQUAD  (N,  F,  X,  T,  MACHEP,  H) 


ARGUMENTS:  N  -  Number  of  variable 

F  -  The  function  to  be  minimized 
X  -  Coordinate  of  the  initial  point 
in  the  quadratic  space 
H  -  Step  size 

COMMON  BLOCKS:  PRCOMM,  Q 


SUBROUTINES  CALLED: 


PRMIN 


CALLED  BY: 


PRAXIS 


NAME: 


SUBROUTINE  PRSORT 


PUR  [’OS  E: 

CALLING  SEQUENCE: 
ARGUMENTS : 

COMMON  BLOCKS: 
SUBROUTINES  CALLED: 
CALLED  BY: 


Sorts  the  elements  of  D(N)  into 
descending  order  and  moves  the 
corresponding  columns  of  V(N,N) 
resulting  in  eigenvalues  and  eigenvecto 

CALL  PRSORT  (M,  N,  D,  V) 

M  -  The  maximum  step  size 
N  -  Number  of  variables 
D  -  Second  difference  array 
V  -  Matrix  of  search  directions 

NONE 

NONE 

PRAXIS 
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NAME: 


SUBROUTINE  QPLOTZ5 


PURPOSE : 

CALLING  SEQUENCE: 
ARGUMENT : 

COMMON  BLOCKS: 
SUBROUTINES  CALLED: 
CALLED  BY: 


Subroutine  QPLOTZ5  calculates  the 
scaling  information  needed  to 
generate  the  format  to  label  the  left- 
hand  side  of  the  printer  plot. 

CALL  QPLOTZ  5  (PDQ) 

PDQ-scaling  factor  for  ordinate  plot 
XRPLOTF ,  XRPLOTQ 
NONE 
PLOT  4 
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NAME: 

PURPOSE: 


CALLING  SEQUENCE: 


FUNCTION  RANDUM 

To  avoid  resolution  valleys,  by  making 
random  steps  on  a  periodic  basis. 

Function  RANDUM  returns  a  random  number  1 

uniformly  distributed  in  (0,1) 

RA  =  RANDUM  (NAUGHT) 


ARGUMENTS : 

COMMON  BLOCKS: 


NAUGHT  -  number  used  to  initialize 
the  random  number  generator 

NONE 


CALLED  BY: 


PRAXIS 


i 
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NAME: 

PURPOSE: 


CALLING  SEQUENCE: 
ARGUMENTS : 

COMMON  BLOCKS: 
SUBROUTINE  CALLED: 


SUBROUTINE  REFLKT 

Subroutine  REFLKT  defines  the  lower 
half  of  the  symmetrical  matrices  T  and 
W,  the  auxiliary  wave  resistance  integrals, 
by  reflection  about  the  diagonal  of  the 
matrix. 

CALL  REFLKT 

NONE 

AUX,  COEFS 
NONE 


CALLED  BY: 


RWAVE 
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NAME: 


SUBROUTINE  RINIT 


PURPOSE: 

CALLING  SEQUENCE: 
ARGUMENT : 

COMMON  BLOCKS: 
SUBROUTINE  CALLED: 
CALLED  BY: 


Subroutine  RINIT  initializes  the  auxiliary 
wave  resistance  matrices,  T  and  W,  to  zero 

CALL  RINIT 

NONE 

AUX,  COEFS 

NONE 

RWAVE 


; 

?  I 
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\ 

! 
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NAME: 


SUBROUTING  RINTEG 


PURPOSE : 


Subrouting  RINTEG  evaluates  the 
integrand  for  T  and  W  functions. 


CALLING  SEQUENCE: 


CALL  RINTEG  (ALFA,  B,  D,  NL0C2,  WTINT , 
SEPCOS,  SQ) 


ARGUMENTS : 


ALFA 

Integrating  variable 

B 

CSTRT1 (I)/XLS 

D 

CSTRT2 (I) /XLS 

NL0C2 

Flag  indicating  the 
presence  of  a  second  strut 

NL0C2 

=  0  if  single  strut 
s  1  if  tandem  strut 

WTINT 

-  Weighting  constant 
for  the  integrand 

SEPCOS 

-  Value  of  the  cosine 
function  in  the 
integrand 

SQ 

-  Value  of  the  integrand 

1 

(02  _  y2) 1/2 


COMMON  BLOCKS: 
SUBROUTINE  CALLED: 
CALLED  BY: 


AUX,  OMEGA,  COEFS,  WORKVL 

BESSJ 

RWAVE 


COMMENTS:  Integrals  used  in  calculations  are  of  the  form: 
J  Yos  +  *  da  r  <*>  da 


>/  a2  - 


f(o) 


/, 


'OS  v  -  '08  YOS  +  *  ■  -  OS 

For  the  first  integral  we  use  the  substitution:  a  =  y 


1  v/  a2  - 


os 


f(a) 


+  C2  where  da  ■  . 


r1  2d£ 

Therefore  the  first  integral  becomes  J  - —  f(Yos  +  C2) . 


2  os 

Where  the  function  f(a)  takes  the  form  of  E  (a)  J  (a)J  (a) 

s  m  n 

where  the  order  of  the  Bessel  functions  is  determined  by  which  auxiliary 
function  is  to  be  evaluated. 


NAME: 


SUBROUTINE  ROUT 


P"RPOSE: 


CALLING  SEQUENCE: 
ARGUMENTS : 


COMMON  BLOCKS: 
SUBROUTINES  CALLED: 
CALLED  BY: 


Calculates  EHP,  and  wave  resistance  coefficients, 
frictional  drag  coefficients,  speed  constants, 
and  prints  the  intermediate  and  final  design 
results 

CALL  ROUT  (TITLE,  EHP,  KEY) 

TITLE  -  label  printed  on  output 

EHP  -  effective  horsepower  of  the  ship 

KEY  -  variable  used  to  skip  or  print  output 

OMEGA,  PHYSCO,  WORKVL 

SUM,  FORMDR 

Function  EHP 


COMMENTS:  The  type  of  output  from  ROUT  is  determined  by 

the  value  of  KEY,  which  is  set  in  the  subroutines 
which  call  EHP.  If  KEY  =  0  (set  in  PEHP)  no 
output  is  printed  and  ROUT  just  does  computations; 
if  KEY  =  1  (set  in  SWATHO)  the  optimization  is 
complete  and  the  final  design  values  are  printed; 
if  KEY  =  2  (set  in  AUXPLOT)  the  intermediate 
design  results  are  printed 


r 
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NAME: 

PURPOSE: 

CALLING  SEQUENCE: 
ARGUMENTS : 


SUBROUTINE  RWAVE 

Subroutine  RWAVE  computes  the  auxiliary  wave 
resistance  functions  T  and  W. 

CALL  RWAVE  (B,  D,  NLOC2) 

B  -CSTRT1 (I) /XLS 

D  - CSTRT2 ( I ) /XLS 

NL0C2  -Flag  indicating  the  presence  of  a 
second  strut 


'{? 


NL0C2-r  if  8in8le  strut 
if  tandem  struts 


COMMON  BLOCKS: 
SUBROUTINES  CALLED: 
CALLED  BY: 


AUX,  OMEGA,  PSI,  COEFS,  PHYSCO,  WORKVL 
RINIT,  SIMPSN,  RINTEG,  BESSJ,  REFLKT 
FUNCTION  EHP 


COMMENTS : 


Integrals  evaluated  by  RWAVE: 


*  (2m-l)  (2n-l)  /. 


da 


mn 


W„ 


^os  a2 


v42  -  Y2 


os  3^(0)  3^(0) 


=  (2m)  (2n)  /  -§■  -  Y,2 


mn 


mn 


os 


os  E  (a)  J  ,  (a)  J_  (a) 

s  2m  2  n 


(2m-l)  (2n-l)  /. 


a2da 


Yos  /a2  -  y2~ 
08 


w„ 


(2m)  (2n)  / 


a2da 


mn 


os 


vS2"  Y 


2 

08 


SB, 


mn 


W, 


SB, 


mn 


(2m-l)  (2n-l)  /“  dot 
(2m)  (2n) 


°s  v/a2  -  y2 
os 


/°°  da 
Y°s 


V8)  J2m-1(8)  J2„-l<6> 


M6>  J2„«> 


V>  VW  J2»-l(0)  Vl<" 


V’  EB<6)  J2„<“>  J2„<S> 


OS 
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NAME: 


SUBROUTINE  SCALE 


PURPOSE: 


GALLIC 3  SEQUENCE: 
ARGUMENTS : 

COMMON  BLOCK: 
SUBROUTINES  CALLED: 
CALLED  BY: 


Creates  current  working  values  X,  by  X  »  S*XI, 
where  S  is  the  scale  factor,  and  XI  is  the 
vector  of  initial  values. 

CALL  SCALE  (S,  N) 

S  -  scaling  value  for  each  of  the  variables 
N  -  number  of  variable 
INITL,  NOCHG,  PHYSCO,  WORKVL 
NONE 

SWATHO,  PEHP,  AUXPLT 
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NAME: 


SUBROUTINE  SIMPSN 


PURPOSE: 

CALLING  SEQUENCE: 


Subroutine  SIMPSN  sets  up  Simpson's  multipliers 
for  numerical  integration  by  Simpson's  rule. 

CALL  SIMPSN  (NPTS ,  SIMP) 


ARGUMENTS : 

NPTS 

SIMP 

COMMON  BLOCKS: 

NONE 

-  Number  of  integration  steps 

-  Array  containing  the  Simpson's  multipliers 


SUBROUTINE  CALLED: 


NONE 


CALLED  BY:  RWAVE 

COMMENT:  These  values  are  used  to  integrate  the  auxiliary 

wave  resistance  integrals,  T  and  W,  for  values 

of  the  integral  between  Y„_  and  Y  _  +  1 

OS  os 


An  example  of  one  of  these  integrals  in  the  range  YQa  to  YQ8+1  is 


+1 


da 

>/a2  -  Y2 
08 


f(a)  ■ 


C  1  ,2dC 

Jo  ^ 2y°a 


+  C 


2 


mos  +  c2) 
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NAME: 


SUBROUTINE  SUM 


PURPOSE:  Subroutine  SUM  computes  the  sums  for  the  wave 

resistance  from  the  Chebychev  coefficients  and 
the  auxiliary  wave  resistance  functions,  T  and  W 


M  N 

Z  I  (A  AT  +  BBW) 

,  ,  m  n  mn  m  n  mn 

m*l  n»l 

COMMENTS:  where  Am  and  Bm  are  the  symmetric  and  antisymmetric  Chebychev 

coefficients  for  the  strut  or  body,  and  Tmn  and  Wm  are  the 
auxiliary  wave  resistance  functions  for  various  hull-strut 
combinations 


CALLING  SEQUENCE: 


CALL  SUM  (SUM1S ,  SUMLB,  SUM1SB,  SUM12S, 
SUM12B,  SUM12SB) 


ARGUMENTS : 


COMMON  BLOCKS: 


SUM1S 

-  Partial 

sum 

for 

strut  1 

SUM1B 

-  Partial 

sum 

for 

body  1 

SUM1SB 

-  Partial 

sum 

for 

Interaction 

between 

strut  1 

and  body  1 

SUM12S 

-  Partial 

sum 

for 

interaction 

between 

strut  1 

and  strut  2 

SUM12SB 

-  Partial 

sum 

for 

interactions 

between 

strut  1 

and  body  2  or 

strut  2 

and 

body  1 

AUX,  COEFS 


SUBROUTINE  CALLED: 


NONE 


CALLED  BY: 


ROUT 


SUMI5 


Figure  4  -  Body-Strut  Configuration 
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NAME: 


SUBROUTINE  VCPRINT 


PURPOSE: 

CALLING  SEQUENCE: 
ARGUMENTS: 


COMMON  BLOCKS: 
SUBROUTINES  CALLED: 
CALLED  BY: 


With  certain  input  options,  this  subroutine 
prints  a  vector  of  given  length  using  two  of 
several  formats 

CALL  VCPRINT  (OPTION,  V,  N) 

OPTION  =  1  Prints  the  Second  Difference  Array 

=  2  Prints  the  Scale  Factors 

=  3  Prints  the  value  of  the  approximating 

quadratic  form 
=  4  Prints  the  value  of  X 

V  =  Matrix  of  search  directions 
N  =  Number  of  variable 

NONE 

NONE 

PRAXIS 


96 


NAME: 


SUBROUTINE  WSURFB 


PURPOSE: 

CALLING  SEQUENCE: 
ARGUMENTS : 

COMMON  BLOCKS: 
SUBROUTINES  CALLED : 
CALLED  BY: 


Determine  surface  area  of  body  of  revolu¬ 
tion  whose  sectional  area  is  given  by 
Chebychev  coefficients 

CALL  WSURFB  (AREA) 

2 

AREA  -  Surface  area  of  a  body  (ft  ) 

COEFS,  PHYSCO,  WORKVL 

FUNCTION  EVAL 

CHEB 
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NAME: 

PURPOSE: 

CALLING  SEQUENCE: 
ARGUMENTS : 

COMMON  BLOCKS: 
SUBROUTINES  CALLED: 
CALLED  BY: 


SUBROUTINE  WSURFS 

Determines  surface  area  of  a  strut  whose  thick¬ 
ness  distribution  is  given  by  Chebychev 
coefficients 

CALL  WSURFS  (AREA) 

2 

AREA  -  Wetted  surface  area  of  a  strut  (ft  ) 
COEFS,  WORKVL 
FUNCTION  EVAL 
CHEB 
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NAME: 


FUNCTION  YINTP 


PURPOSE: 

CALLING  SEQUENCE: 


COMMON  BLOCKS: 


Function  YINTP  interpolates  through  a  set  of 
discrete  data,  using  quadratic  interpolation. 

YTP  =  YINTP  (XZ,  X,  Y,  N) 

XA  -  Point  to  be  interpolated 
X  -  Array  of  ordinate  data 
Y  -  Array  of  abscissa  data 
N  -  Number  of  data  points 

NONE 


SUBROUTINE  CALLED: 


NONE 


CALLED  BY: 


FORMDR 


COMMENTS : 


Uses  linear  interpolation 


I. 
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APPENDIX  A 

OBJECTIVE  AND  PENALTY  FUNCTIONS 
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APPENDIX  A  -  OBJECTIVE  AND  PENALTY  FUNCTIONS 

The  program  SWATHO  is  designed  to  determine  the  geometry  of  a  SWATH  ship  which 
yields  the  minimum  effective  power  for  a  given  speed.  Thus  the  objective  function 
is  the  effective  power.  However,  this  is  an  over  simplification  of  the  problem, 
in  that  the  solution  to  the  problem  must  also  satisfy  a  number  of  constraints.  That 
this  is  so  is  intuitively  obvious  if  one  considers  that  the  ship  with  minimum  power 
for  a  given  speed  is  that  ship  with  zero  length  and  volume,  and  therefore  zero 
resistance.  Thus  an  obvious  constraint  is  that  the  ship  contain  at  least  some 
minimum  volume. 

The  subroutine  PRAXIS,  which  is  used  to  perform  the  minimization,  is  a  program 
designed  for  use  in  unconstrained  minimization  problems.  This  dictates  that  a 
modified  objective  function  be  employed,  which  takes  into  account  the  fact  that 
constraints  are  or  are  not  violated.  To  meet  this  requirement,  an  external  con¬ 
straint  method  was  employed,  where  the  functions  g.(x)  were  defined  to  be  greater 

t  li  ^  t  Vi 

than  zero  when  the  i  n  constraint  was  violated,  and  less  than  zero  when  the  i 

constraint  was  satisfied.  Using  these  g^'s,  the  modified  objective  function  was 

defined  as: 

PEHP  -  EHP(1  +  Max(g  (x) , 0) ) . 
i 

In  all  of  the  above  definitions,  the  vector  jc  is  the  vector  of  parameters  defining 
the  ship's  geometry. 

This  definition  of  the  objective  function  results  in  the  addition  of  no  penalty 
to  the  objective  function  internal  to  the  region  where  the  constraints  are  satisfied, 
and  results  in  the  addition  of  a  penalty  exterior  to  the  region  where  the  constraints 
are  satisfied.  This  gives  rise  to  the  description  of  this  as  an  external  constraint 
method. 

There  are  a  total  of  nineteen  penalty  functions  which  are  applied  to  the  pro¬ 
blem  of  minimizing  the  resistance  of  a  SWATH  ship.  These  constraints/penalty  func¬ 
tions  are  as  follows: 

1.  Minimum  displacement  constraint 
gl  -  (DISPMN-DISP)a1 

where  DISPMN  is  the  minimum  ship  displacement  as  input  by  the  user  and 
DISP  is  the  current  ship  displacement  at  that  iteration  of  the  calculations. 
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2.  Minimum  draft-diameter  ratio  constraint 

g2  *  [Draft  Factor  x  BDIA  -  (HS  +  BDIA)J  c*2 

where  the  Draft  Factor  is  the  minimum  draft  to  diameter  ratio,  BDIA 
is  the  body  diameter  at  half  of  the  length  of  the  body  and  HS  is  the  draft 
of  the  strut. 

3.  Maximum  draft  constraint 

g3  =  (HS  +  BDIA  -  DRFTMAX)a3 

where  DRFTMAX  is  the  maximum  draft  of  the  ship  as  input  by  the  user. 

4.  Maximum  beam  constraint 
g4  *  (BS  +  BDIA  -  WMAX)«4 

where  BS  is  the  separation  of  the  hull  centerlines  and  UMAX  is  the 
maximum  beam  of  the  ship  as  specified  by  the  user. 

5.  Hull  overlap  constraint 
g5  -  (BDIA  -  BS)a5 

6.  Strut-hull  Nose  clearance  constraint 
g6  -  (XFWD  +  CSTRT  +  XLS/2  -  XLB/2)a6 

where  XFWD  is  the  minimum  distance  from  the  body  leading  edge  to  the 
strut  leading  edge,  CSTRT  is  the  distance  of  the  strut  forward  of  the 
body  centerline,  XLS  is  the  length  of  the  strut  and  XLB  is  the  length  of 
the  body. 

7.  Strut -hull  tail  clearance  constraint 
g?  -  (XAFT  -  XLB/ 2  +  XLS/2  -  CSTRT) a? 

where  XAFT  is  the  minimum  distance  from  body  trailing  edge  to  strut 
trailing  edge. 

8.  Minimum  transverse  GM  constraint 
g8  -  [miNGMT  -  (KB  +  BMT  -  KG)]ag 

where  MINGMT  is  the  minimum  transverse  GM,  KB  is  the  distance  of  the 
center  of  buoyancy  from  the  baseline,  BMT  is  the  transverse  metacenter  and 


103 


KG  is  the  distance  of  the  center  of  gravity  from  the  baseline. 

9.  Minimum  longitudinal  GM  contraint 
g9  -  [mINGML  -  (KB  +  BML  -  KG)]ot9 

where  MINGML  is  the  minimum  longitudinal  GM  and  BML  is  the  longitu¬ 
dinal  metacenter. 

10.  Body  offsets  greater  than  zero  contraint 

g10  “  /_!  dx[lb<x>l  "  b(x)]“l0 

where  b(x)  represents  the  body  offsets 

11.  Strut  offsets  greater  than  zero  constraint 
gu  -Ji  dx[|t(x)|  -  t(x)]an 

where  t(x)  represents  the  strut  offsets 

12.  Maximum  length-diameter  ratio  check 
g12  -  (XLB  -  LODMAX  x  BDIA)a12 

where  LODMAX  is  the  maximum  body  length  over  diameter  ratio. 

13.  Minimum  length-diameter  ratio  check 
g13  -  (LODMIN  x  BDIA  -  XLB)a13 

where  LODMIN  is  the  minimum  body  length  over  diameter  ratio. 

14.  Strut  minimum  thickness  constraint 
g14  -  (TSMIN  -  TSMAX)a14 

where  TSMIN  is  the  minimum  strut  thickness  at  50%  chord  and  TSMAX  is 
the  strut  thickness  in  the  current  iteration  of  the  calculations. 

15.  Strut  thickness-body  diameter  constraint 
g15  -  (TSMAX  -  BDIA)ai5 

16.  Positive  strut  length  contraint 
g16  -  (-XLS)a16 

17.  Positive  body  length  constraint 

m  (-XLBjoijy 
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18.  Maximum  body  offset  constraint 

g  a  *  (b  -  1.2  BDIA)a. R 
18  max  18 

where  b  is  the' maximum  body  offset.  This  constraint  insures  a 
max 

reasonable  value  for  the  maximum  body  offset. 

19.  Maximum  strut  thickness  constraint 

*19  -  (t.ax  -  ‘-2  TS«AX)o,9 

where  t  is  the  maximum  possible  value  of  strut  thickness.  This 
max 

constraint  insures  a  reasonable  value  for  the  maximum  strut  thickness. 

In  all  of  the  above  penalty  functions,  the  constraints  a  are  chosen  so  that 

1  4 

a  constraint  violation  of  ten  percent  results  In  a  value  of  around  10  . 

Of  the  above  penalty  functions,  eleven  are  constraints  Imposed  by  the  program 
user  to  insure  that  the  optimum  vehicle  meets  some  minimum  set  of  criteria  such  as 
minimum  displacement  and  transverse  and  longitudinal  GM.  These  user  specified 
constraints  are  the  constraints  numbered  1,  2,  3,  4,  6,  7,  8,  9,  12,  13,  and  14. 

The  other  eight  constraints  are  constraints  Imposed  by  the  program  to  insure 
that  the  geometry  of  the  vehicle  is  reasonable.  Examples  of  these  would  be  con¬ 
straints  to  insure  that  the  lengths  of  the  body  and  strut  are  positive,  and  con¬ 
straints  to  insure  that  the  offsets  of  the  strut  and  body  are  all  positive.  The 
numbers  of  these  constraints  are  5,  10,  11,  15,  16,  17,  18,  and  19. 
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APPENDIX  B  -  COMPUTATIONAL  PROCEDURE  FOR  SWATH  RESISTANCE  PREDICTION 


The  main  purpose  of  the  SWATHO  program  is  to  minimize  the  power  required  to 
propel  a  SWATH  ship  in  a  calm  sea  at  constant  speed.  For  a  ship  moving  at  speed  V 
and  experiencing  a  total  drag  force  RT  the  effective  horsepower  required  is: 

P  -V 

E  550  ‘ 

The  total  drag  is  comprised  three  components:  frictional  resistance  (Rp) , 
wave-making  resistance  (R^i,  v*rm  drag  (R^) .  That  is: 

T  RF  +  ^  +  RFM* 

The  frictional  resistance  and  :'rm  drag  are  caused  by  the  motion  of  the  hull  through 
a  viscous  fluid.  The  wave-making  resistance  is  due  to  the  energy  that  must  be 
supplied  by  the  ship  to  the  wave  system  created  on  the  free  surface.  However,  wave¬ 
making  resistance  and  form  drag  are  usually  grouped  together  under  the  title  of 
residuary  resistance. 

FRICTIONAL  RESISTANCE 

Frictional  resistance  is  the  single  largest  component  of  the  resistance  of  a 
SWATH  ship.  Frictional  resistance  is  calculated  in  the  traditional  fashion  of  naval 
architects,  using  the  ITTC  model  ship  correlation  line  to  determine  a  frictional 
resistance  coefficient  (Cp) .  However,  for  a  SWATH  ship,  the  frictional  drag  of  the 
hulls  and  struts  are  calculated  separately  based  on  their  respective  Reynolds 
numbers.  Thus  we  have  that  the  frictional  resistance  coefficient  is  given  by  the 
formula: 


CF  " 


0.075 


{log1Q  Rn  -  2) 2 

and  the  frictional  resistance  of  the  hulls  and  struts  are  given  by: 

RF  =  ^  SHull  VZ  CF 
Hull  Hull 


and 


^  SStrut  VZ  CF  * 
Strut  Strut 
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In  calculating  the  total  frictional  resistance,  a  correlation  allowance  (CA)  of 
0.0005  is  also  included  in  the  total,  so  we  have  that: 


RF 


+  Rr- 


Hull 


+  HP  (S. 


Strut 


Hull 


SStrut)V2 


FORM  DRAG 

Form  drag  (R^)  is  usually  defined  as  the  viscous  component  of  the  drag  due  to 
the  shape  of  body,  i.e.,  the  difference  between  the  total  viscous  drag,  and  the 
viscous  drag  of  the  equivalent  flat  plate  of  the  same  length.  However,  in  the  case 
of  the  SWATH  resistance  programs,  the  form  drag  has  been  determined  in  a  much  more 
empirical  fashion. 

SWATH  form  drag  has  been  defined  as  the  difference  between  the  experimentally 
determined  residuary  resistance  of  a  hull  form  and  the  theoretically  derived  wave 
resistance  for  that  same  hull  form.  Such  a  difference  was  calculated  for  both 
SWATH  3  and  SWATH  4  based  on  bare  hull  resistance  results, "*■  and  plotted  as  a  func¬ 
tion  of  strut  speed-length  ratio.  A  curve  was  then  faired  through  the  data. 

Figure  3  (p.  69).  Due  to  the  scatter  in  the  data  a  decision  was  made  not  to  allow 
the  form  drag  coefficient  go  below  0.0005;  this  can  be  seen  on  Figure  3.  The  form 
drag  coefficient,  determined  from  Figure  3,  is  converted  to  full  scale  form  drag, 
using  the  traditional  formula  and  the  total  wetted  surface: 

■'em  '  *>  SV’  CFM  ' 


wave-making  resistance 

The  wave-making  resistance  of  a  ship  is  not  easily  predicted  from  empirical 
relations  or  from  gross  ship  features.  Experience  has  shown  that  two  ships  of 
similar  form  may  differ  significantly  in  their  measured  total  resistance  due  to 
differences  in  the  wave-making  resistance  component.  Wave-making  resistance  is  a 
function  of  the  fine  details  of  ship  form,  and  is  best  studied  through  mathematical 
analysis . 

Lin  and  Day^  investigated  the  problem  of  wave-making  resistance  for  twin-hull 
ships.  The  mathematical  problem  is  formulated  within  the  context  of  the  linearized 
thln-ship  theory.  In  their  investigations,  the  potential  flow  around  such  a  ship 
is  represented  by  sheet  distributions  of  sources. 

The  computer  program  SWATHO  documented  in  this  user's  manual  is  devoted 
primarily  to  the  implementation  of  numerical  procedures  for  the  prediction  of  the 
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wave-making  resistance  of  SWATH  ships  following  the  theory  developed  by  Lin  and 
Day.  The  remainder  of  this  appendix  will  be  devoted  to  the  discussion  of  the 
essential  details  involved  in  the  numerical  procedure  utilized  in  the  wave-making 
resistance  predictions. 

Based  on  the  analysis  of  Lin  and  Day,  the  wave  resistance  is  expressed  as 
follows : 

RS  ■  (t  <*TVos) 

M  M 

x  ^  S  {ASm  ASn  TSmn  +  BSm  BSn  WSmnl , 
m=l n=l  1  f 


RB 


R 


SB 


„  2  v  M  M 

2TrpgAn_,_  , 

LsYos  /  (  ABm  ABn  TBmn  +  BBm  BBn 

m=l  n=l 


(27rpgApT) 


M  M 

zz 

m=l  n=l 


Sm  Bn  SBmn 


+  B 


D 

Sm  Bn 


w. 


SBmn 


W 

SBmn 


} 

} 


where  T,  Lg,  and  hg  are  the  maximum  thickness,  length,  and  draft  of  the  strut, 
respectively;  A0,  Lg,  and  hfi  are  the  maximum  section  area,  length,  and  depth  of 
submergence  of  the  axis  of  the  body. 

In  the  above  equations,  R  ,  R  ,  and  R  represent  the  wave  resistance  due  to 

b  d  bfi 

one  strut,  one  main  body,  and  the  intersection  between  strut  and  main  body,  respec¬ 
tively.  Hence,  the  wave  resistance  of  a  SWATH  ship  becomes 


■Hi  -  2(rs  +  rb  +  rsb>- 


Also  in  the  above  equations,  Agm  and  Bgm  are  the  Chebychev  coefficients  for  the 

strut,  A  and  B  are  the  Chebychev  coefficients  for  the  hull,  and  T_„  ,  T_„  , 

Bm  Bin  SMmn  SBmn 

^Bmn’  WSmn’  WSBmn’  anc*  WBmn  are  t^e  aux^^^ary  wave  resistance  functions.  The 

auxiliary  wave  resistance  functions  are  defined  as  follows: 


T 

Smn 

(2m-l^ (2n-l) 

Wc 
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(2m) (2n) 
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jJ2n)(a)  J2n(ot) 
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and 


^  ,  ,  2b  2  1-5 — 

D  -  1  +  cos  —  — 

S  Tos 


2 

os 


Es  =  1  -  e 


gLS 

fos  "  (2U7)’ 
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(2ra-l)(2n-l) 
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where  J  (a)  is  a  Bessel  function  of  the  first  kind,  defined  as  follows: 
m 

2  r 11,2 

J  (a)  =  —  I  cos(a  sin  t)  cos(mt)  dt. 
m  tt  I 

0 


NUMERICAL  INTEGRATION  METHODS  FOR  AUXILIARY  WAVE 
RESISTANCE  FUNCTION  EVALUATION 

\ 

The  range  of  integration  for  the  auxiliary  wave  resistance  functions  is  divided 
into  three  intervals  of  integration:  lYos»  YQS  +  1],  a  "central  region"  and  the 
"tail"  region.  A  detailed  discussion  of  these  intervals  and  the  integration  proce¬ 
dure  for  each  interval  follows. 

INTEGRATION  PROCEDURE  FOR  THE  INTERVAL  [y  ,  y  +1] 

'os’  'os 

In  the  region  of  yog  the  integrands  do  not  behave  well,  but  they  are  neverthe¬ 
less  integrable.  This  is  due  to  the  presence  of  ./a2-  y2  in  the  denominator.  This 

v  'os 

may  be  easily  handled  by  utilizing  the  substitution: 

a  =  y  +  S2 

os 

where 

da  =  2£d£. 

This  substitution  will  be  utilized  in  the  region  [y  ,  y  +11.  It  eliminates  the 

°  'os  'os 

square-root  singularity  by  giving 


Investigating  the  behavior  of  the  integrands  as  functions  of  X,  in  the  region 
0  <  t,  <  1  shows  that  they  are  very  well  behaved  and  a  simple  numerical  integration 
scheme  with  about  30  points  will  yield  excellent  results. 


Ill 


INTEGRATION  PROCEDURE  FOR  THE  CENTRAL  REGION 


The  central  region  of  the  numerical  integration  starts  at  a  =  y  +  1  and  con¬ 
tinues  until  Eg(6)  reaches  a  very  small  fraction  of  its  initial  value. 


EB<6i„U>  ‘  EB 


<Y„s+1)YOB 


OS 


—  c 

The  fraction  is  given  as  10  ,  so 


2  hfi 

Bmax2  ( 

exp 

s 

Y0B  j 

( 

,  hB 

exp  | 

'  lb 

Yo.\b  ( 

Hence, 


The  justification  of  this  choice  for  a  ,  and  the  effect  of  the  choice  of  E,  will 

max 

be  discussed  in  the  next  section. 

In  order  to  implement  an  effective  and  economical  numerical  integration  scheme, 
it  is  necessary  to  investigate  the  behavior  of  the  integrands  within  this  central 
region.  The  requirement  is  to  provide  at  least  several  points  of  the  integrand  for 
each  full  cycle  of  its  oscillation.  Hence,  a  conservative  estimate  of  behavior  will 
be  an  estimate  of  more  rapid  oscillation  than  actually  exists. 

The  Bessel  functions  are  approximated  as 

v>«4¥ 

for  |z|»l  and  jz|»|v|. 
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for  smaller  Z,  the  behavior  will  be  less  oscillatory,  so  this  represents  a  con¬ 
servative  estimate  for  all  Z.  This  assumption  results  in 


f  4  +  1] 


J  (a)Jv(a)  = 


even 


[cos(2a+'J')  -  l]  ^  odd 


where  the  case  of  (y-v)  odd  will  never  occur  due  to  the  form  of  the  integrands. 
The  other  combinations  of  Bessel  functions  are  found  as 


J  (a)J  (B) 
m  n 


tt/SB 


cos(a+B) 


Jm(8)J„(W  -  -ig-  cos(26) , 


where  the  phase  angles,  a  slow  oscillation  in  the  first  case,  and  a  constant  term 
in  the  second  case,  have  all  been  ignored.  Since  3  is  generally  greater  than  a, 
but  of  the  same  order  of  magnitude,  a  conservative  estimate  gives  the  variation  of 
all  of  the  above  as 


J  J 


m  n 


4  -  4- 


The  term  [ 4x2-y^g  1  behaves  as  [a]  for  a  »  YQS»  and  is  better  behaved  for  smaller 
a.  This  conservative  estimate  will  be  used  for  the  occurrence  of  this  term  at  the 
beginning  of  each  integrand.  However,  the  substitution  will  not  be  made  in  the  term 


D(a)  =  cos 


2b 


—  -L-a/a^-r 

'os  S 


T 

OS 


which  occurs  in  all  of  the  two-hull  integrands,  since  that  would  result  in  consid¬ 
erably  less  economical  numerical  integration. 

The  terms  Eg(B)  and  Eg (a)  contain  decreasing  exponentials  and  are  not  oscil¬ 
latory  and  are  quite  well-behaved,  and  so  will  not  affect  the  spacing  of  points  for 
the  numerical  integration. 
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The  single-hull  terms  (TSmn>  WSmn,  TBnm,  WRmn,  TSBmn,  WSBmn)  are  similar  but 
I  2  2b  i 

without  the  cos  - —  — — ajn2--v2  term. 

LYos  lS  v  'os  J 

Based  on  the  behavior  shown  above,  the  most  rapid  oscillation  of  the  integrands 
is  like 


cos 

In  the  vicinity  of  a  *  a. 


cos(f(a)]  behaves  like  cos 


a=a 


and  hence 
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If  the  accuracy  desired  from  the  numerical  integration  scheme  requires  P  points  per 
cycle,  then  the  step  size  S  (distance  between  points)  in  the  vicinity  of  a  =  a  is 


S  = 


2tt 

P 


2b  2S2-y*s 
LS  ^2-Yos  - 


This  step  size  will  decrease  with  increasing  a  but  will  be  smaller  for  larger  hull 

centerline  spacing  b.  Hence,  if  computations  are  made  for  various  values  of  b  for 

the  same  sample  points,  the  step  size  should  be  based  on  the  largest  value  of  b. 

The  numerical  integration  scheme  chosen  for  the  central  region  is  a  modified 

Simpson's  rule  in  which  all  steps  are  not  equal.  For  this  method,  a  step  size  is 

computed  at  a  given  a  (starting  first  with  Y0g  +  1)*  using  the  last  equation.  Then 

an  integral  is  taken  over  an  interval  equal  to  two  of  these  step  sizes,  using  a 

three-point  Simpson's  rule.  A  new  step  size  is  computed  for  the  a  at  the  end  of 

that  interval,  and  the  process  is  repeated  until  a  is  reached  or  exceeded.  The 

nidx 

concept  may  be  visualized  as  shown  below. 


integrind 


SurfSOfB  Simpscre 

ftul*  Rule 


This  economizes  on  samples  (points)  where  the  integrand  does  not  vary  quickly,  and 
uses  more  densely  spaced  points  where  necessary. 

All  of  the  approximations  made  in  this  section  have  been  only  for  the  purpose 
of  choosing  the  step  sizes  for  the  numerical  integration.  The  samples  taken,  which 
are  multiplied  by  appropriate  weights  added  to  obtain  the  numerical  evaluation  of 
the  integral,  are  all  found  using  the  actual  integrands. 
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INTEGRATION  PROCEDURE  FOR  THE  TAIL 


It  should  be  noted  that  the  criterion  for  the  end  of  the  "central  region,"  and 
hence  the  beginning  of  the  "tail,"  is  that  Eg ($)  reach  a  very  small  percentage  of 
its  original  value.  This  will  insure  that  any  integrand  which  includes  E  (8)  will 
have  long  ceased  to  contribute  significantly  to  the  result.  The  body-alone  and 
strut-body  integrals  (Tg,  Wg,  T^,  W12g,  Tgg,  Wgg,  T^g  and  W12SB)  can  therefore  be 
assumed  complete  at  the  conclusion  of  the  "central  region"  numerical  integration. 

The  strut-alone  terms  (Tg,  Wg,  T^g.  and  W^2g)  are  not  complete,  however.  The 
integrand  for  these  terms  decreases  only  like  a-4.  The  utilization  of  the  "central 
region"  numerical  integration  scheme  until  this  integrand  were  very  small  would  be 
quite  accurate  but  exceptionally  expensive,  especially  due  to  the  constantly  de¬ 
creasing  (with  increasing  a)  step  size  required  for  that  scheme.  Hence,  a  more 
economical  scheme  must  be  developed,  accepting  the  resulting  loss  of  accuracy. 

The  single-hull  strut-alone  terms  (Tg  and  Wg)  are  easily  handled  since  they  do 
not  include  the  highly  oscillatory  term 


D(a) 


cos 


os 


2b 

L_ 


a/ot2-Y 


2 

OS 


The  behavior  of  the  integrands  for  these  terms  is  like 


Eg(ot)  cos(2a) , 


tra 


which  allows  the  use  of  a  simple  numerical  integration  scheme  with  constant  step 
size 

s  -  y 


where  P  is  the  desired  number  of  points  per  cycle. 

The  above  method  is  still  not  very  economical  due  to  the  large  region  over 
which  the  integral  must  be  taken,  due  to  the  relatively  slowly  decreasing  a-4 
envelope  in  which  the  integrand  oscillates.  A  more  economical  method  is  to  approx¬ 
imate  the  integrand  asymptotically,  and  analytically  integrate  the  approximate 
integrand  from  a  given  value  to  infinity.  The  obvious  simplification  is  to  assume 
that  Eg(a)  *  1  ,  which  is  an  excellent  approximation,  becuase  Eg(ot)  approaches  unity 
as  Eg(8)  approaches  zero.  The  approximation  is  hence  assured  by  the  condition  that 
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determined  the  end  of  the  "central  region."  In  addition  to  Eg(a),  the  Bessel 
functions  must  also  be  approximated.  We  use  the  approximation 


Jv(a) 


f  2  /„  VTT  TT  \ 

=  yjm-  C0s\z  "  T  -  ~) 


for  |a|  »  1  and  |a|  »  |v| . 

If  |a|  is  not  enough  greater  than  |v|  at  the  beginning  of  the  tail,  a  special 
"base  of  tail"  numerical  integration  is  done  for  Tg  and  Wg  using  step  size 
S  =  u/P  (as  discussed  on  the  previous  page) .  This  is  done  over  a  region  from  the 
start  of  the  "tail"  until  a  is  large  enough  for  the  above  asymptotic  approximation 
to  yield  an  accurate  result.  The  above  asymptotic  approximation  gives 


Jy(v)Jv(y)  =: 


-  [cos(2a+'{')  +  1]  -x-  ever 

ira  2 

~  [cos(2a+4*)  -  1]  odd 

ira  2 


where  the  case  of  (y-v)  odd  will  never  occur  due  to  the  form  of  the  integrands. 

The  integrands  therefore  consist  of  a  constant  term  and  an  oscillatory  term.  The 
contribution  of  the  constant  term  to  the  integral,  using  the  above  approximations. 


arc  sinl 


a.  ffci:’./Vv2~v2 


a  ^os 


m- 


v/a2-Y2 
__A  *  os 

27rYosaA2  ’ 


When  is  much  greater  than  YQg»  this  expression  will  not  provide  an  accurate 

solution  since  it  will  be  the  difference  of  two  very  similar  values.  For 

R  =  y  /a.  »  1,  the  result  is 
os  A 


— ir-  R  4-  4-  R3  +  .4;  R5  +  .  . 

2<s  6  40 


1  +  ^  R2  +  .  .  .  . 

3tto  *  10  j 
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The  simplified  form 


1 

3™A 

will  therefore  be  accurate  to  about  2n  decimal  places  if  R  is  10  n.  Likewise,  the 
exact  solution  will  lose  2n  decimal  places  in  the  subtraction  of  the  two  parts  when 
R  is  10  n.  Hence,  for  a  computer  with  about  16  decimal  place  accuracy,  the  simpli¬ 
fied  form  is  usfed  whenever  R  <  10-4. 

The  analytic  integration  of  part  of  the  integrands  for  Tg  and  Wg  shown  in  the 
preceding  paragraph  depends  on  a  large  enough  for  sufficiently  accurate  approx¬ 
imation  of  the  Bessel  functions.  Hence,  the  simple  numerical  integration  shown 
must  be  utilized  to  fill  in  the  gap  if  the  "central  region"  ends  before  such  a 
sufficiently  large  aA  is  reached. 

The  oscillatory  components  of  the  integrands  for  Tg  and  Wg  in  the  region  of 
analytic  integration  of  the  asymptotic  approximations  have  been  ignored  in  this 
formulation,  as  their  contribution  is  assumed  to  be  negligible. 

The  integrals  for  and  W^2g  were  not  carried  past  the  end  of  the  "central 

region"  due  to  their  highly  oscillatory  integrands,  which  made  the  contribution 
from  further  computation  both  uneconomical,  and  not  critical  to  the  final  values  of 
T12s  and  W12g.  This  is,  however,  not  as  good  an  approximation  as  those  previously 
noted  in  the  above  development.  As  a  result,  errors  of  as  much  as  2%  can  be 
expected  in  the  computation  of  wave  drag  due  to  interference  between  the  two  struts. 


APPENDIX  C 

RELATIONSHIPS  BETWEEN  CHEBYCHEV  SERIES  AND  SWATH 
HULLFORM  COEFFICIENTS 


APPENDIX  C  -  RELATIONSHIPS  BETWEEN  CHEBYCHEV  SERIES  AND  SWATH 

HULLFORM  COEFFICIENTS 


CHEBYCHEV  SERIES 

In  the  wave  resistance  integrals,  geometry  of  strut  and  body  are  represented 
by  a  special  form  of  Chebychev  series.  In  SWATHO,  the  Chebychev  coefficients  are 
■pproximated  from  various  statistical  moments  of  strut  and  body.  The  Chebychev 
series  representation  of  hull  and  body  geometries  are  discussed  below.  Let 

x  =  sin  6,  --j-  <_  6  <_  —  . 

Define  the  fundamental  functions  of  the  Chebychev  series  as  follows: 

U  (x)  *  cos(2m-l)0 
m 

=  cos[(2m-l)sin-1x] 


V  (x)  =  sin  2m8 
m 

=  sin{2m  sin-1x],  m  =  1,  2,  .  .  .  ,  m. 


The  strut  half-thickness  and  the  body  sectional-area  functions  can  be  repre¬ 
sented  by  finite  sums  of  the  fundamental  Chebychev  series  as  follows: 


and 


M 


t(x)  [A 

m=l 


U  (x) 
sm  m 


B  V  (x)l 
sm  m  J 


M 

=V  [a  cos(2m-l)8  +  B  sin  2m9l 
L  sm  sm  J 

m=l 


M 


A(x)  [A. 

=  1 


bmVX> 


Bu  V 
bm  m 


(x)] 


m= 


M 

“S  [AbmCOS(2,n~l)0  +  Bbmsin  2m03 

m=l 
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Through  use  of  the  orthogonality  of  the  series,  the  following  inversion 
formulas  are  obtained: 


Sm 


.if1 

"  J- 1 

-  -  f 

71  J 

•fir 


dx 
tt/2 


t(x)U(x) 

m 


d0t (sin  0)  cos[(2m-l)0]. 


t(x)V  (x) 
m 


Sm 


■  Hi  -  ^ 

■ 1  /l-x 


•tt/2 


2  r1 

—  I  d0t (sin  0)  sin(2m0) , 
-tt/2 


and 


ABm 


2_  f1,/2 

77  J  d0  A(sin  0)  cos(2m-l)0, 
-tt/2 


BBm 


d0  A(sin0  )  sin  2m0. 


—  tt/2 


INVERSE  PROCEDURE  FOR  DETERMINING  THE  CHEBYCHEV  COEFFICIENTS 

The  geometric  coefficients  of  a  SWATH  ship's  form  can  be  shown  to  be  closely 
related  to  the  Chebychev  coefficients  of  the  strut  and  body. 

These  coefficients  of  form  are  easily  obtained  at  the  preliminary  design  stage 
and,  thus,  can  be  used  to  approximate  the  Chebychev  coefficients.  In  the  following 
example,  the  waterplane  coefficient,  C^p,  can  be  used  to  find  AC1  . 

Wr  SX 

By  definition,  the  waterplane  area  is 

Ls/2 

Ay  =  2  /  dx  t (x) 

J-Ls/2 
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where  t(x)  is  the  thickness  function  of  the  strut.  Let 


Ls  Ls  1 

x  =  ~2~  ^ ■  dx  =  ~2~  d^»  C  =  iTyj  x*  =  tmax  * 

S 


Thus , 


c 


t  /  dCt(C) 
2  .  max  * 


2  max 


/_!  «[£)  V, 


(C>  +  B  V  (?) 
m  sm  m 


1 


-T  W  |  £  //  «„«>«  +  £  Bsm  /  ‘  vm(Odfl. 


Making  the  substitution. 


x  *  sin  0  ,  dx  =  cos  0  d0, 


we  obtain: 


Ls  [  it,  r  ff/2 

^  2  tmax  |  Ef  Asm  J_^,2 


d0  cos(2m-l)0cos  0 


M  TT/2 

+  EBsm  / 

1  ■'-tt/2 


d0  sin  2m0  cos  0 


The  second  integral  is  odd  and,  therefore,  it  is  identically  zero. 

M  /.  tt/2 


L  t 

,  s  max 
\ - 2 - 


n  r 

Ea  f 

sm  J  . 


t  M  r 

s  max  .  / 

~~2  4^AsmJ- 


— tt/2 
tt/2 


d0{cos(2m-l)  0  cos  0  ] 


2 -  4^"smJ_Tr/2  d0 [cos(2m-2)0  +  cos  2m0] 


Hence 
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1 


L  t 
s  max 

2 


L  t 
s  max 

2 


A 

sm 


A 

sm 


1 

2m-2 


sin(2m-2)6  + 


2m 


sin  2m0 


Tt/2 

0 


sin(m-l)iT  sin  nm 
2(m-l)  2m 


These  integrals  are  identically  zero  except  for  m=l  where  the  first  term  contri¬ 
butes  tt/2. 


AW  = 


L  t 
s  max 

2 


7T 

2 


A 


si- 


The  waterplane  coefficient  is  defined  as 


C 


WP 


L  t 
s  max 


Therefore 


4C 


SI 


WP 
7 r 


This  derivation  shows  how  the  waterplane  coefficient  is  related  to  the  Chebychev 
coefficient,  Ag^.  Similarly  we  can  prove, 

4C„ 


B1 


where  Cp  is  the  body  prismatic  coefficient, 


In  the  SWATHO  program  the  maximum  number  of  terras  of  the  Chebychev  series  is 
three.  The  Chebychev  coefficients  are  determined  by  the  following  formulas  for  the 
strut : 


4C 


WP 


SI 


S2 


AS1(1  "  16  CIW)’ 
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AS3  1  "  AS1  AS2’ 


bsi  =  4clcf  '  Asr 


BS2  ■  °' 


BS3  “  °' 


and  the  body  coefficients  are  determined  as  follows: 


abi  =  T  V 


AB2  1  "  AB1’ 


AB3  =  °' 


BB1  “  4CLCF  *  AB1! 


BB2  °> 


where 


BB3  ”  °‘ 


C  =  Waterplane  area  coefficient 

W  r 


CWP  =  '  V'  where  \  is  waterplane  area  of  one  strut. 

CT „„  =  Waterplane  moment  coefficient 
LCr  ^ 

CLCF  =  *  LS^ 

C-r„  =  Waterplane  inertia  coefficient 


"IW 


CIW  =  V(AW  ’  V 


Cp  =  Body  prismatic  coefficient 


JI,CB 


C„  =  V_/(Av  •  L_)  where  VD  is  displaced  volume  of  one  body. 
P  B  X  d  B 

=  Body  moment  coefficient 

MXB/(VB  *  LB^ 


All  moments  and  the  moment  of  inertia  are  taken  about  the  mid-length  of  the 
respective  strut  or  body. 


12A 


DTNSRDC  ISSUES  THREE  TYPES  OF  REPORTS 

1.  DTNSRDC  REPORTS.  A  FORMAL  SERIES,  CONTAIN  INFORMATION  OF  PERMANENT  TECH 
NICAL  VALUE.  THEY  CARRY  A  CONSECUTIVE  NUMERICAL  IDENTIFICATION  REGARDLESS  OF 
THEIR  CLASSIFICATION  OR  THE  ORIGINATING  DEPARTMENT. 

2.  DEPARTMENTAL  REPORTS,  A  SEMIFORMAL  SERIES,  CONTAIN  INFORMATION  OF  A  PRELIM 
INARY,  TEMPORARY,  OR  PROPRIETARY  NATURE  OR  OF  LIMITED  INTEREST  OR  SIGNIFICANCE 
THEY  CARRY  A  DEPARTMENTAL  ALPHANUMERICAL  IDENTIFICATION. 

3.  TECHNICAL  MEMORANDA.  AN  INFORMAL  SERIES,  CONTAIN  TECHNICAL  DOCUMENTATION 
OF  LIMITED  USE  AND  INTEREST.  THEY  ARE  PRIMARILY  WORKING  PAPERS  INTENDED  FOR  IN 
TERNAL  USE.  THEY  CARRY  AN  IDENTIFYING  NUMBER  WHICH  INDICATES  THEIR  TYPE  AND  THE 
NUMERICAL  CODE  OF  THE  ORIGINATING  DEPARTMENT.  ANY  DISTRIBUTION  OUTSIDE  DTNSRDC 
MUST  BE  APPROVED  BY  THE  HEAD  OF  THE  ORIGINATING  DEPARTMENT  ON  A  CASE  BY  CASE 
BASIS. 


